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PART   I. 
THE  Uoji  OF  THE  NONLINEAR  PROGRAMMING  TECHNIQUE  POP-II  TO 
OPTIMIZE  a  MULTI-EFFECT  MULTI-ST^SE   (S4EMS) 
SE&WATEH  DISTILLATION  PLANT 


1.0   INTRODUCTION 

In  this  part  of  the  report  a  Multi-Effect  Multi-Stage  (MEMS) 
seawater  distillation  plant  was  optimized  by  applying  the  nonlinear 
programming  technique  known  as  POP-II.   Several  concepts  basic  to  the 
understanding  of  nonlinear  programming  will  be  presented  first. 
1.1   Some  Basic  Concepts  and  Definitions  Useful  in  Nonlinear  Programming 

Unfortunately  no  single  algorithm  exists,  with  the  exception  of  a 
systematic  exhaustive  search,  for  solving  the  general  nonlinear  program- 
ming problem  which  has  an  objective  function  of  the  form 

S  =  f(x1,  x2 ,  x   )   m   f(x)  (la) 

X  —  (.Xj  ,  x~,  ....,  X  ) 

and  constraints  of  the  form 


,00 


t>-  (lb) 


i  =  1,  2,  ....,  m 

It  should  be  apparent  that  there  is  no  single  method  which  is  best 
for  attacking  the  above  problem.   Many  times  the  best  method  for  attacking 
a  problem  will  depend  a  great  deal  upon  the  exact  functional  forms  of  the 
objective  function  and  the  constraints.   For  example,  if  all  the  constraints 
are  equality  constraints,  the  methods  given  by  Hadley  (1)  can  be  applied. 

A  few  special  cases  of  the  general  nonlinear  programming  problem  will 
be  given  in  the  sections  that  follow.   Fortunately,  most  problems  of 
practical  interest  can  be  handled  by  one  or  more  of  the  special  cases  or 
by  more  general  techniques  such  as  POP-II. 


. 


Before  proceeding  further,  certain  definitions  which  are  important  to 
nonlinear  programming  will  be  introduced. 

A  bounded  region  of  n  dimensional  space,  sometimes  called  n-space  or 
hyperspace,  is  a  portion  of  hyperspace  which  has  been  closed  off  or  sur- 
rounded.  Hypersurfaces  are  used  to  form  the  bounded  region.   The  concept 
of  a  bounded  region  is  best  visualized  in  terms  of  2  or  3  dimensional 
space.   Figure  la  shows  a  portion  of  the  first  quadrant  which  is  bounded 
by  plane  curves.   Figure  lb  shows  a  portion  of  the  first  octant  which  is 
bounded  by  a  spherical  surface. 

A  point  which  is  inside  a  bounded  region  is  called  an  interior  point. 
A  point  on  the  boundary  is  a  surface  point.   A  point  outside  of  the 
bounded  region  is  called  an  exterior  point.   Interior  and  surface  points 
are  feasible  points  whereas  an  exterior  point  is  said  to  be  infeasible. 

The  quadratic  form  is  useful  in  nonlinear  programming  in  that  by 
examining  the  quadratic  form  of  a  function  it  is  possible  to  make  important 
statements  about  the  nature  of  the  function.   A  function  is  said  to  be 
strictly  concave  if  its  quadratic  form  is  negative  definite  over  all  of 
Euclidean  space.   For  a  function  with  one  variable  this  corresponds  to  the 
case  in  which  there  is  a  unimodal  maximum.   On  the  other  hand,  if  the 
quadratic  form  of  the  function  is  positive  definite,  the  function  is  said 
to  be  strictly  convex.   This  corresponds  to  a  function  with  a  unimodal 
minimum  in  the  one  variable  case.   Convex  and  concave  functions  of  a  single 
variable  are  shown  in  Figures  2a  and  2b,  respectively. 

Recall  that  if  the  quadratic  form  of  a  function,  h(x),  is  positive 
definite,  then  the  quadratic  form  of  -h(x)  is  negative  definite.   This 
implies  that  if  h(x)  is  convex,  then  -h(x)  is  concave  and  conversely.   Or, 
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if  h(x)  is  unimodal  with  a  maximum,  then  -h(x)  is  unimodal  with  a  minimum. 
That  is,  it  is  possible  to  change  a  minimization  problem  into  a  maximiza- 
tion problem  or  a  maximization  problem  into  a  minimization  problem.  There- 
fore, one  may  talk  exclusively  about  either,  maximization  or  minimization 
processes. 

The  following  relations  also  hold,  where 

(Convex  Function)    f%    |  1x^(1-  X  )x2  j  £  \f  ^x,  )+(l-  A  )  f^)     (2a) 

(Concave  Function)  i%    ("Ax^I-  X  )x2  ]  >    XfjfcjHU-X)  i^x^  (2b) 

One  way  of  thinking  about  the  significance  of  relation,  Equation  (2a),  is 
to  consider  the  problem  of  approximating  the  convex  function  pictured  in 
Figure  2a  by  a  straight  line  through  any  two  points  on  the  curve.   Observe 
that  the  straight  line  lies  above  the  function  it  is  approximating.   That 
is,  it  over  estimates  the  function.   Similarly,  the  straight  line  approxi- 
mation to  a  concave  function  always  under  estimates  the  function  it  is 
approximating.   Relations  given  by  Equations  (2a)  and  (2b)  can  be  genera- 
lized to  higher  spaces  by  replacing  x  by  the  s-dimensional  vector,  x  - 
(x^>  x9'  -•••'  xs)*   *n  t'le  general  case  of  s-space  these  relations  are 
orten  used  to  define  convex  and  concave  functions. 

The  important  use  of  the  definitions  of  convexity  and  concavity  is  to 
determine  whether  or  not  an  extremum  is  a  local  or  a  global  extremum.   For 
the  case  of  minimization  of  an  objective  function,  a  local  minimum  which 
occurs  within  a  closed  convex  constraint  set  is  a  global  minimum  if  the 
objective  function  is  convex.   A  local  maximum  is  a  global  maximum  if  the 
objective  function  is  concave  and  the  constraints  are  closed  and  convex. 
Note  in  both  cases  the  constraint  set  must  be  closed  and  convex. 
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For  many  problems,  especially  those  with  large  numbers  of  constraint 
equations,  it  is  quite  a  laborious  task  to  determine  if  the  objective 
function  is  convex  or  concave  and  if  the  constraint  equations  are  convex. 
Therefore,  what  is  generally  done  is  to  solve  the  problem  without  examining 
the  object  function  and  constraints  with  respect  to  convexity  and  con- 
cavity.  Simulation  and/or  a  search  technique  are  then  used  to  verify  the 
results.   If  this  reveals  a  solution  point  which  gives  a  better  extremum 
than  the  previous  solution,  the  new  point  is  used  as  a  starting  point  for 
applying  the  particular  nonlinear  programming  method  again.   The  danger  of 
a  local  extremum  should  be  kept  well  in  mind,  but  for  many  problems  a 
local  optimum  is  better  than  none  at  all. 

The  careful  reader  will  also  note  that  it  is  possible  for  an  objective 
function  to  be  convex  over  parts  of  the  constraint  region  and  concave  over 
other  parts.   This  case  is  shown  in  Figure  3. 

Furthermore,  some  functions,  such  as,  the  one  given  by  equation  (3), 
may  be  neither  convex  nor  concave. 

S  =   C    a.x.  (3) 

i-1   1  X 

1.2  The  Kuhn-Tucker  Conditions 

The  Kuhn-Tucker  conditions  are  useful  because  they  give  some  insight 
into  nonlinear  programming  theory  and  they  enable  one  to  determine  whether 
or  not  a  point  which  has  been  found  by  a  computational  procedure  is  an 
optimum  point. 

Kuhn  and  Tucker  (2)  were  the  first  to  derive  the  conditions  which  bear 
their  names.  Carr  and  Howe  (3)  also  give  a  very  readable  derivation  of  the 
conditions. 
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The  derivation  of  the  Kuhti-Tucker  conditions  given  here  is  due  to 
Wilde  and  Beightler  (4). 

Consider  the  problem: 

min:     S  =  f(x)  (4a) 

subject  to  the  constraints 

%i    (x)  i.  b.  (4b) 

There  is  no  loss  in  generality  by  using  constraint  equation  (4b),  in 
place  of  constraint  equation  (lb)  since  constraint  equation  (4b)  can  be 
obtained  from  constraint  equation  (lb)  by  introducing  more  slack  variables 
and  constraints. 

Constraint  equation  (4b)  can  be  charged  to  an  equality  constraint  by 
introducing  the  slack  variable,  u,  which  is  squared  to  insure  that  the 
squared  quantity  remains  positive. 

g.  (x)  -  b.  -  u.2  -  0  (5) 

A  Lagrangian  function  can  be  formed 

L  =  f(x)  -   £  *i  J  g,  (x)  -  b.  -  u.2  (6) 


Then,  three  of  the  four  Kuhn-Tucker  conditions  which  are  necessary 
conditions  for  a  minimum  are: 

in 

J  l  -  0  ■  Jf(x)  -   t_.   Ai  J  gi(x),  j  =  1,  2 s   (7) 


cJx. 

h 


=  0 


1,  2,  ...,  m       (S) 


g,-(x) 
cUi 
,)l    =  0  =  2/L  a.,    i  =  1,  2,  ...,  m  (9) 

dui 
The  fourth  necessary  Kuhn-Tucker  condition  for  a  minimization  problem  as 

given  by  Wilde  and  Beightler  (4)  requires  that  the  Lagrange  Multipliers  be 


non-negative, 

^  i  -  0,   I  -  1,  2,  ....,  m  (10) 

If  a  point  does  satisfy  the  Kuhn-Tucker  conditions,  it  may  or  may  not 
be  a  minimum.   However,  if  it  does  not  satisfy  them,  it  can  not  possibly 
be  a  minimum.  If  both  the  objective  function  and  the  constraints  are 
convex,  the  Kuhn-Tucker  conditions  are  sufficient  conditions  for  a  global 
minimum.  The  Kuhn-Tucker  conditions  in  themselves  do  not  provide  a  com- 
putational procedure  for  finding  optimal  solutions  to  nonlinear  programming 
problems. 

2.0  NONLINEAR  PROGRAMMING  METHODS 

2.1  The  Classical  Calculus 

For  some  nonlinear  programming  problems,  classical  calculus  techniques 
can  be  used  to  find  the  optimum  extreme  point.  The  steps  for  finding  the 
optimum  extreme  point  given  the  objective  function,  equation  (la),  and  the 
constraints,  equation  (lb),  are  as  follows:  Ignoring  the  inequality  con- 
straints, find  the  stationary  points  of  the  objective  function.  Determine 
which  if  any  of  the  stationary  points  are  within  the  feasible  region.  For 
those  points  in  the  feasible  region  determine  at  which  point  the  objective 
function  takes  on  an  extreme  value.  Search  the  boundaries  of  the  feasible 
region  to  determine  if  the  interior  extreme  point  is  a  global  extremum. 

An  example  illustrating  this  procedure  is: 

Find  the  maximum  value  of  the  objective  function 

2 
S  ■  IOx-l  +  20x2  +  x.x.  -  2xL  -  2x2  (11) 

which  is  subject  to  the  constraints 

0  -  xx  i   7  (12) 

0  -  x2  -  8  (13) 
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The  stationary  points  and  the  values  of  the  objective  function  at  these 
points  are  given  in  Table  1.   For  this  problem  the  maximum  is  at  the  in- 
terior point.   Note  however,  that  if  the  goal  would  have  been  to  minimize 
the  objective  function,  the  optimal  point  would  have  been  on  the  boundary. 

Several  mathematical  programming  algorithms  are  extentions  of  the 
simplex  method  of  linear  programming.  The  reluctance  to  abandon  the  sim- 
plex method  may  be  partially  explained  by  noting  that  large  sums  of  money 
have  been  invested  in  the  development  of  computer  codes  to  solve  large 
scale  linear  programming  problems  since  1947  when  the  simplex  method 
started  its  rapid  development.  Furthermore,  experienced  people  are  avail- 
able who  can  cast  a  nonlinear  programming  problem  into  a  linear  program- 
ing format. 

Two  classes  of  nonlinear  programming  problems  which  can  be  solved  by 
modified  simplex  algorithms  are  separable  and  quadratic  programming  problems. 

2.2  Separable  Programming  Problems 

The  separable  programming  technique  can  be  applied  when  the  constraints 

are  linear  and  it  is  possible  to  separate  the  objective  function  into  a  sum 

of  functions  each  of  which  is  a  function  of  only  one  independent  variable. 

The  general  separable  programming  problem  has  the  following  form. 

Objective  function: 

f(x)       f     t.    (x  )  (14) 

i=l 

Constraint  equations: 

t    aij  xj  (=)    V  i  =  l m  <15) 

j=l 

2.3  Quadratic  Programming   Problems 

When  the  objective  function  has  the  quadratic  form 


n 


Table  1.    Stationary  Points  and  Values  of  the  Objective 
Function  for  the  Classical  Calculus  Example 


Location 

of 

Stationar 

y 

Point 

X-i 

x2 

5 

interior 

4.0 

6.0 

80.0 

Boundary 

0.0 

5.0 

50.0 

Boundary 

7.0 

6.8 

63.1 

Boundary 

4.5 

8.0 

72.5 

Boundary 

2.5 

0.0 

12.5 

12 


f(x)  -  E    a.x   +    T      HC   xixi  (16) 

1-1   x  x       1-1     j-1   «  l  J 


and  the  constraints  are  linear, 

£    a.,  x.   U  /   b.;  i  =  1 ,  m  (17) 

j-i   1J  J  UJ    l 

quadratic  programming  methods  can  be  used  to  extremize  the  objective 
function. 

Separable  programming,  quadratic  programming,  and  methods  for  solving 
other  types  of  nonlinear  programming  problems  are  discussed  in  detail  by 
several  authors  (1,  3,  4).   In  this  work,  little  emphasis  will  be  placed 
on  methods  for  solving  problems  which  have  specific  objective  and  constraint 
equation  forms  because  more  general  algorithms  exist  for  solving  more  general 
types  of  mathematical  programming  problems. 

2.4  The  Process  Optimization  Program  (POP-II) 

For  many  of  the  nonlinear  programming  algorithms,  no  "general",  easy  to 
use,  computer  program  (computer  code)  exists  for  implementing  the  calculations 
for  large  scale  problems.   Because  this  is  a  serious  disadvantage  to  using 
these  methods,  there  has  been  a  great  impetus  to  develop  a  general  computer 
program  to  handle  a  wide  variety  of  nonlinear  optimization  problems,  which 
would  require  a  minimum  of  effort,  knowledge,  and  time  on  the  users  part. 
One  computer  code  for  doing  this  which  has  been  quite  successful  is  the 
Process  Optimization  Program  II  (POP-II)  which  has  been  developed  by  Smith 
(5),  of  IBM's  System  Research  Institute. 

POP-II  uses  a  truncated  Taylor  series  to  obtain  a  sectionally  linear- 
ized linear  programming  problem  from  the  nonlinear  programming  problem. 
To  start  the  computations,  the  processes'  performance  and  objective  equations 
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are  linearized.   A  linear  programming  problem  is  formed  from  these 
linear  equations.   The  solution  of  the  linear  programming  problem,  hope- 
fully, gives  values  of  the  independent  variables  that  are  closer  to  the 
optimum  values.   The  above  procedure  is  repeated  until  the  optimum  values 
of  the  independent  variables  have  been  determined. 

POP-II  is  a  nonlinear  programming  technique  which  is  easy  to  use. 
All  that  is  necessary  for  using  the  technique  is  a  mathematical  model  of 
the  process  or  system  which  the  user  wishes  to  optimize  and  a  mathematical 
statement  of  the  optimization  objective  (objective  function).   The  model 
must  be  constructed  so  that  the  processes'  dependent  variables  can  be 
calculated  once  the  independent  variables  have  been  specified.  If  the 
user  can  provide  a  model  which  will  calculate  the  dependent  variables  given 
the  independent  variables,  then  POP  may  well  be  a  good  technique  to  use. 

Some  of  the  terms  discussed  in  this  section  may  be  easier  to  understand 
if  one  refers  to  Appendix  II  where  the  computer  output  for  a  POP  problem 
is  given. 

The  process  model  or  simulation  program  which  the  user  provides  is 
named  SUBROUTINE  MODEL.   SUBROUTINE  MODEL  must  be  programmed  in  such  a 
way  that  each  of  the  processes'  dependent  variables,  the  y.  values,  is 
written  as  a  function  of  one  or  more  of  the  processes'  independent  variables, 
the  x  values,  and/or  a  function  of  one  or  more  of  the  previously  defined 
dependent  variables.   In  effect  each  dependent  variable  must  be  written 
as  a  function  of  one  or  more  of  the  independent  variables. 

The  sectionally  linearized  linear  programming  technique  used  by  POP 
is  described  in  the  following  paragraphs. 

Consider  the  objective  function  and  the  equality  constraints  of  the 
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generalized  nonlinear  programming  problem  given  by  equations  (la)  and 
(lb). 

S  »  f(xy,   x, ,  Xg)  (18) 

gp  (*l.  *2 ,  xs)  =  bp,  p  =  1,  2,  ...,  n  (19) 

If  the  inequality  constraints  are  ignored  for  the  present,  the  number  of 
independent  variables  for  the  nonlinear  programming  problem  posed  by 
equations  (18)  and  (19)  is  equal  to  (s-n).   Let  the  number  of  independent 
variables  be  r,  that  is,  r  is  equal  to  (s-n). 

Let  the  i  th  dependent  variable  be  given  by  ?..  However,  y1,  is 
reserved  for  the  objective  function  S.  It  is  now  possible  to  rewrite 
equations  (18)  and  (19)  in  the  form 

f(y,,  y3,  .....  yn+1,  xv   *2,  *r)  (2°) 


and 


where 


gp  (y2.  y3,  •••..  yn+1.  V  x2>  Kr)  "  V 

p  -  1,  2 ,  n  (21> 


Equations  (20)  and  (21)  represent  respectively,  the  general  functional 
forms  of  the  objective  function  and  the  equality  constraints.   In  most 
instances,  an  individual  g  will  not  be  a  function  of  all  of  the  y^, 
i  =  2,  3,  n+1  and  x  ,  y  =  1,  2,  ....,  r.   For  purposes  of  program- 
ming a  SUBROUTINE  MODEL,  each  equality  constraint  is  solved  for  a  particular 
y.  value.   The  equality  constraints,  having  been  solved  for  a  particular 
y.  value,  are  next  arranged  in  an  orderly  fashion,  such  that  each  y.  is  a 
function  of  the  other  dependent  variables  y2,  y^,  ....,  Y^_i  an&  the 
independent  variables,  that  is, 

y2  =  h2  (x) 

y  ■  h3  (h2(x),  x)  =  h^(x) 
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y.   =  h.  (h,<x),  h3(x) ,  h1_1  (x),  x)   =  h.  (x)      (22) 


=  h  ,,  (h.  (x),  h,(x),  ,  h  (x),  x)  -  h  . ,  (x) 


'n+1   "n+1  v"2 


n+Iv 


where  the  abbreviated  notation  x  '  (x.,  x,m  ....  x  )  is  used.   After  the 

values  of  y,,  y,,  ....,  y  ,,  have  been  determined,  y.  is  calculated  using 

equation  (20). 

The  dependent  variables  given  by  equation  (22)  may  be  approximated  in 

a  small  region  near  a  point  x  ,  denoted  by  x   =   (x   ,  x   ,  ....,  x  ),  by 

o    o         o 
a  truncated  Taylor  Series  if  a  coordinate  system  is  defined  such  that  the 

distance  from  x  in  the  i  th  direction  is  denoted  by  S  x. .   Thus  we  can 


,       j=r         o>   hi 
y.  (xq+  Sx)  ~  h.(x)      +  C    Sx    ^x 


(23) 


lj  <-,  ....5  n 


Equation  (23)  is  linear  in  Sx.   and  can  be  simplified  to  give 

~        j=r         r 
y.  (x  +   x)  -   e.  +   T3    a-'    Jx.  , 
j  =  l 
where 

i  =  2,  3 ,  n+1 


(24) 


e.=  h.  (x) 


(25) 


and 


lj 


(26) 


A  similar  approximation  holds  for  y.  in  the  small  region  near  x  . 
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The  set  of  linear  equations  obtained  by  linearizing  the  nonlinear 
equations  may  be  used  to  form  a  linear  programming  (LP)  problem.   The 
LP  problem  is  solved  in  a  small  region  near  x  •   In  most  instances  the 
solution  will  give  a  new  x  point  which  gives  an  improved  value  of  the 
objective  function. 

When  inequality  constraints  are  imposed  they  are  introduced  into 
the  LP  problem  by  the  POP  program  after  the  user  has  inserted  the  maximum 
and/or  minimum  limits  of  a  variable  on  the  input  data  cards.   The  POP 
program  handles  the  inequalities  by  automatically  introducing  the  neces- 
sary slack  variables  which  are  required. 

The  size  of  the  small  region  near  x  in  which  the  LP  problem  is 
solved  may  be  controlled  by  user-specified  constants  which  are  called 
move  limits.  The  j  th  move  limit  is  the  maximum  step  size  that  is 
permitted  for   5x.. 

Since  the  y.  values  obtained  in  the  linear  programming  solution  may 

deviate  from  the  true  y.  values,  the  new  x  values  are  used  in  SUBROUTINE 
'i       '  o 

MODEL  to  calculate  the  true  y.  values.   The  user  must  also  supply  error 
limits  for  the  dependent  variables.   These  error  limits  are  used  by  POP 

to  minimize  the  deviations  that  occur  from  the  true  y.  values. 

J  x 

The  procedure  described  up  until  now  is  called  a  loop  and  is  repeated 
until  an  extreme  value  of  the  objective  function  is  found.   In  the  event 
that  an  extreme  value  is  not  found,  the  procedure  is  terminated  when  the 
maximum  number  of  optimization  loops  which  the  user  specifies  has  been 
exceeded. 

Many  nonlinear  optimization  techniques  require  the  calculation  of 
derivatives  of  the  same  type.   A  significant  feature  of  POP  is  that  the 
values  of  the  first  partial  derivatives  required,  the  a. .  values,  are 
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calculated  automatically  by  a  subroutine  of  POP.   The  central  difference 

technique  which  is  used  to  calculate  the  partial  derivatives  is  given  as 

h. 


U 


<)XJ 


£X1     |  jc  I  AXi      (27) 
2  Ax. 


where  Che  Ax.   values  are  specified  by  the  user.   This  technique  is 
illustrated  graphically  in  Figure  4  for  the  case  where  h.  is  a  function 
of  one  independent  variable,  x, .  Note  that  h.  is  shown  as  a  smooth  con- 
tinuous curve  in  Figure  4.   Obviously,  if  the  function  whose  first  partial 
derivative  is  to  be  computed  by  a  numerical  method  has  one  or  more  discon- 
tinuities for  example,  then  the  calculated  derivative  may  be  a  very  poor 
approximation.   This  is  a  disadvantage  inherent  in  the  use  of  POP.   How- 
ever, on  the  other  hand,  techniques  which  require  analytical  expressions' 
for  first  partial  derivatives  require  logical  statements  in  their  computer 
programs  for  handling  piecewise  smooth  functions  (see  Appendix  I  of  Part  II 
of  this  work  for  an  example.) 

As  an  optional  feature,  the  subroutine  which  calculates  the  partial 
derivatives  can  also  calculate  new  move  limit  values  as  the  optimization 
progresses.   The  process  is  called  making  an  adaptive  move  limit  calculation. 
This  feature  or  operating  mode  of  POP  may  be  used  with  some  SUBROUTINE 
MODEL'S  to  reach  the  optimum  in  a  fewer  number  of  loops.   However,  with 
some  SUBROUTINE  MODEL'S  the  adaptive  move  limit  calculation  reduces  the 
values  of  the  move  limits  so  much  that  POP  shuts  off  before  reaching  the 
optimum.   Smith  (5)  gives  a  detailed  description  of  the  adaptive  move 
limit  operating  mode. 

In  writing  the  SUBROUTINE  MODEL  used  to  simulate  the  MEMS  plant, 
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ax, 


Slope    of    £S 


h,(x) 


h. 


'lx,  +AX, 
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x^-Ax, 


Fig.  4.   Central  difference -approximation  to  the  slope  of 
the  function  h.  (x.)  at  the  point  x.  . 
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several  procedures  were  found  to  help  speed  the  programming  and  debugging 
processes.   Since  it  is  thought  that  these  procedures  will  be  helpful  in 
general  programming,  they  will  be  mentioned  here. 

In  SUBROUTINE  MODEL  the  j  th  independent  variable  is  denoted  as  X(J). 
The  i  th  dependent  variable  is  denoted  as  Y(l).  The  i  th  constant  is 
denoted  as  CONST(I).   If  large  numbers  of  the  symbols  X(J),  Y(I),  and 
CONST(I)  are  used  in  a  SUBROUTINE  MODEL,  it  may  be  difficult  to  remember' 
what  each  symbol  represents.  Consequently,  it  may  be  easier  to  debug 
SUBROUTINE  MODEL  by  writing  it  in  terms  of  easily  recogni2able  symbols. 
The  independent  variables  which  are  written  in  terms  of  X(J)  variables 
at  the  beginning  of  the  program.   At  the  end  of  SUBROUTINE  MODEL  the  Y(I) 
variable  can  be  defined  in  terms  of  the  easily  recognizable  symbols  which 
are  used  in  order  to  make  the  programming  easier. 

From  experience  with  the  MEMS  plant  problem,  it  was  also  determined 
that  it  is  not  necessary  to  define  all  dependent  variables  as  Y(l)  vari- 
ables.  For  example,  the  horsepower  requirements  in  the  MEMS  plant  problem 
were  not  defined  by  Y(I)  variables.  Since  the  program,  without  modification, 
has  a  maximum  limit  of  50  Y(I)  variables,  this  procedure  may  also  help 
spread  out  the  Y(I)  variables  and  thus  increase  the  problem  size  that  can 
be  handled. 

The  programming  procedure  described  above  will  require  more  computer 
cards.  However,  from  knowledge  of  the  problem  studied,  it  is  felt  that  the 
additional  computer  cards  in  SUBROUTINE  MODEL  will  not  have  a  significant 
effect  on  the  time  required  for  execution  of  the  POP  program.   Furthermore, 
the  time  from  problem  inception  to  completion  may  be  greatly  decreased. 

One  difficulty  inherent  in  the  use  of  POP,  which  also  occurs  in  the 
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use  of  many  other  techniques,  is  the  possibility  of  finding  a  relative  or 
local  extreraum  rather  than  a  global  extremum.   If  there  are  several  local 
extreme  points  in  the  feasible  answer  space,  POP  will  probably  find  the 
one  which  is  nearest  the  starting  point.  This  difficulty  has  been  called 
"nearsightedness"  by  Baumol  (6,  3).   The  way  that  "nearsightedness"  is 
usually  overcome  is  to  start  the  program  at  several  widely  spaced  feasible 
points  and  then  observe  what  happens.   Another  problem  which  may  be  en- 
countered when  SUBROUTINE  MODEL  is  large  and/or  has  complex  performance 
equations  is  that  it  may  be  difficult  to  find  a  feasible  starting  point. 
However,  to  aid  the  user,  POP  will  double  the  move  limits  on  successive 
loops  until  a  feasible  point  has  been  found  if  one  can  be  found  at  all. 
A  diagram  giving  the  steps  required  to  solve  a  nonlinear  programming 
problem  using  POP  is  given  in  Figure  4a. 

One  of  the  important  features  of  POP  is'  the  matrix  which  is  printed 
out  after  POP  has  reached  an  optimum  point  or  optionally  after  each  loop. 
Each  element  of  this  matrix,  which  is  called  the  DYDX  MATRIX,  is  one  of 
the  a.,  values  given  by  equation  (27).   Specifically,  the  element  a.,  is 
the  number  used  in  the  i  th  row  of  the  LP  tableau  with  the  variable  x;. 
For  the  i  th  row  the  y  variable  is  y.. 

Kith  the  DYDX  MATRIX,  it  is  a  simple  matter  to  perform  an  incremental 
sensitivity  analysis  to  determine  where  to  concentrate  effort  to  improve 
the  process  being  optimized.   In  addition,  inspection  of  the  first  row  of 
the  matrix,  that  is,  the  row  for  the  objective  variable  y1 ,  will  give  a 
good  indication  of  the  nature  of  the  response  surface  near  the  x  point 
for  which  the  matrix  holds.   This  is  especially  true  if  the  A  x.  values 
used  are  small. 
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Derive  process  performance 
equations  and  establish 
objective  function. 

Equations  (20)  and  (21) 


Identify 
dependent 
variables  : 

and 
independent 
variables 


x.  >  X(J) 


Arrange  performance  equations 
and  write  as  PORTRAH  state- 
ments to  obtain 

SUBROUTIHE  KODEL 


Specify  PCP  input  data  such  as 
move  limits,  dependent 
variable  error  limits,  and 
Initial  values  of  X(J) * 


Insert  SUBROUTINE  KODSL  and 
input  data  into  the  POP  deck 


?.un   x^up  starting  with  initial 
X(j)  values 


Calculate  partial  derivatives: 


Form  and  solve  LP  problem  for 
new  X(j)  values  and  Y(l) 
values  which  may  be  incorrect 
due  to  linearization  errors 


Insert  new  X(j)  values  into 
SUBROUTINE  MODEL  to  calculate 
orrect  Y(l)  values 


Fig.  'ta.   The  steps  required  to  solve  a  nonlinear  programming 
problem  using  POP-II. 
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Each  coefficient  of  the  DYDX  MATRIX  can  be  examined  with  respect  to 
sign  and  size.   Consider  the  case  when  the  coefficient  for  the  i  th  row 
and  the  j  th  column  is  -a.  This  indicates  that  if  the  independent  variable 
of  the  j  th  column  is  increased  by  one  unit  then  the  dependent  variable 
of  the  i  th  row  will  be  decreased  by  approximately  a_   of  its  units.   This 
analysis  holds  true  only  for  a  small  region  around  the  x  point  for  which 
the  matrix  holds.  Furthermore,  the  above  statements  are  not  strictly  true 
for  variables  that  appear  in  equality  constraints. 

If  each  element,  a.  .,  of  the  objective  variable  row  of  the  matrix  is 
a  very  small  number,  then  small  charges  in  the  independent  variables  will 
not  change  the  objective  function  significantly.  If  this  situation  occurs 
at  the  optimum,  this  is  called  a  flat  optimum  and  the  response  surface  is 
said  to  be  well  behaved  near  the  optimum. 

However,  if  one  of  the  elements,  a    0f  tne  objective  variable  row 
is  very  large,  then  only  a  small  change  in  the  independent  variable  of 
that  element,  x  ,  will  cause  a  large  change  in  the  objective  function,, 
This  is  the  most  interesting  type  of  problem  as  it  may  be  a  constraint  which 
causes  this  situation  to  happen.   In  the  event  that  it  is  a  constraint  on  a 
variable,  which  causes  the  element  of  the  objective  variable  row  to  be 
large,  one  should  look  into  the  possibilities  for  removing  or  relaxing  the 
constraint  to  some  extent. 

As  an  example,  consider  a  problem  in  which  reactor  temperature  appears 
as  an  independent  variable.   Furthermore,  assume  an  upper  limit  is  placed 
on  the  reactor  temperature  because  the  reactor  cannot  withstand  high  tem- 
peratures.  Now  if  the  system  is  optimized  by  POP  and  it  is  found  that  the 
reactor  temperature  element  of  the  y.  row  of  the  DYDX  MATRIX  is  the  only 
large  element  in  the  row,  and  furthermore  if  this  is  due  to  the  upper 
temperature  imposed,  then  obviously  searching  for  construction  materials 


23 


which  can  withstand  higher  temperatures  should  be  seriously  considered. 
Admittedly  this  attack  will  necessitate  resolving  the  problem  with  new 
cost  data. 

3.0   OPTIMIZATION  OF  A  MEMS  PLANT  USING  POP-II 

An  example  problem  of  a  Multi-Effect  Multi-Stage  (MEMS)  seawater 
distillation  plant  which  demonstrates  the  use  of  POP-II  in  process  design 
will  be  given  next. 
3.1  Process  Description 

Fan,  et.  al.  (7)  give  an  excellent  description  of  this  process  which 
is  included  here  for  convenience. 

Figure  5  illustrates  a  three-effect  multistage  flash  system.   In 
order  to  facilitate  the  discussion,  some  critical  locations  in  the  system 
are  denoted  by  letters,  Z,  A,  B,  C,  D,  etc.  and  the  system  is  divided  into 
various  sections  which  are  denoted  by  HR-1,  R-l,  HR-2,  etc.   The  first 
effect  consists  of  a  brine  heater  H-l,  a  heat  recovery  section,  HR-1,  and 
a  cooling  section,  R-l.   The  second  effect  consists  of  a  brine  heater  H-2, 
a  heat  recovery  section,  HR-2,  and  a  cooling  section,  R-2.   The  third 
effect  consists  of  a  brine  heater  H-3,  a  heat  recovery  section,  KR-3,  and 
a  cooling  section,  R-3.  The  section  between  locations  B  and  C  serves  a 
double  purpose;  it  is  the  cooling  section  for  the  first  effect,  R-l,  and 
the  brine  heater  for  the  second  effect  H-2.   Similarly,  the  section  between 
locations  E  and  F  serves  as  the  cooling  section  for  the  second  effect,  R-2, 
and  the  brine  heater  for  the  third  effect,  H-3.   Sea  water  is  used  as  a 
coolant  in  R-3. 

F  and  L  represent  the  flow  rate  of  feed  brine  and  flashing  brine, 
respectively*   The  feed  brine  and  recycle  brine  are  referred  to  together  as 
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the  non-flashing  brine  stream.  T,,  T.  and  T  represent  respectively  the 

r   j      c 

temperature  of  the  flashing  brine,  non-flashing  brine,  and  condensate, 
respectively.   Subscript  notation  will  be  used  to  indicate  the  location. 
For  example,  (T-)„,  (T.)„  and  (T  )„  represent  the  temperature  of  the 
flashing  brine  at  location  F,  temperature  of  the  non-flashing  brine  at 
location  B  and  the  temperature  of  the  condensate  at  location  H  respectively. 
R  ,  R,  and  R,  represent  the  recycle  flow  rate  in  the  first,  second  and 
third  effect,  respectively,  and  W.,  W,  and  W,  represent  the  condensate 
produced  in  the  first,  second,  and  third  effects,  respectively.  R,  repre- 
sent the  cooling  water  (sea  water)  used  in  the  third  effect. 

The  sea  water  feed  is  heated  in  R-3  and  then  acidified  and  degasified 
to  remove  C0_  and  other  dissolved  gases.  After  being  heated  successively 
in  KR-3,  H-3,  HR-2,  H-2,  HR-1  it  is  mixed  with  recycle  brine  R,  to  form  a 
brine  stream  which  is  heated  in  brine  heater  H-l  and  then  introduced  into 

the  first  effect  as  the  flashing  brine  (L)  . 

A 

The  flashing  brine  at  location  C  is  divided  into  two  streams.  One 
stream,  (L)„,  is  fed  into  the  second  effect  and  the  other  stream,  R, ,  is 
recirculated  by  a  recycle  pump,  J.,  heated  in  HR-1  and  then  mixed  with  the 
feed  stream  at  the  mixing  point  M-,  As  has  been  described,  the. combined 
stream  is  heated  in  brine  heater  H-l  and  introduced  into  the  first  effect 
as  flashing  brine  (L),.   Similarly  the  flashing  brine  at  locations  F  and  K 
are  divided  into  (L)„  and  R„ ,  and  (L)„  and  R„ ,  respectively.   (L)   is  fed 
into  the  third  effect,  and  (L)   is  discharged  from  the  third  effect  as  the 
reject  brine  from  the  system.   Streams  R„  and  R„  are  recirculated  by  pumps 
J,  and  J,  respectively,  heated  in  sections  HR-2  and  H-2,  and  HR-3  and  H-3, 
respectively,  mixed  with  (L)„  and  (L)   at  mixing  points  M,  and  M,  respective- 
ly, and  introduced  to  the  second  effect  and  the  third  effect  as  (L)B  and 
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(L)  ,  respectively. 

The  feed  brine  and  the  recycle  brine  are  heated  in  each  stage  by  the 
water  vapor  evaporated  from  the  flashing  brine  in  that  stage.   It  is  possible 
to  arrange  the  flow-system  so  that  the  temperatures  of  the  feed  brine  and  the 
recycle  brine  are  equal  at  any  location.   In  the  following  discussion,  such 
an  arrangement  is  assumed.   As  has  been  described,  the  feed  brine  and  the 
recycle  brine  are  referred  to  jointly  as  the  non-flashing  brine  and  its 
temperature  is  denoted  by  T..   The  recycle  brine,  R  ,  which  is  a  part  of 
the  flashing  brine  at  location  C,  is  introduced  into  the  condensing  chamber 
at  location  B  where  it  becomes  a  part  of  the  non-flashing  brine  stream. 
Therefore,  the  following  relation  should  hold. 

<Vb"(Vc  <28> 

Similarly,  we  can  write 

(T.)E  =  (Tf)r  (29) 

and 

(Vh=(Vk'    .  (30) 

Mixing  is  thermodynamically  irreversible  when  two  solutions,  which 
differ  in  temperature  and/or  composition,  are  mixed  together.   The  solutions 
mixed  at  the  mixing  points,  M^,  M,  and  M  ,  have  differences  in  composition; 
however,  by  suitably  locating  point  B,  the  temperature  of  the  recycle 
solution  R,,  (T.)fi,  can  be  adjusted  to  (Tf)  ,  the  temperature  of  stream  ' 
(L)c,  and  in  this  way  the  thermodynamic  irreversibility  due  to  mixing  can 
be  minimized.   In  this  study  isothermal  mixing  is  at  M.  ,  M  ,  and  M  and  the 
heat  of  mixing  of  the  sodium  chloride-water  system  is  neglected.   Thus, 
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(Vb=  (Vc=  (TfV  (3« 

Similarly, 

(T.)E-  (Tf)F=  (Tf)G.  (32) 

The  unit  enthalpy  of  dilute  salt  solutions  is  assumed  to  be  a  function 
of  temperature  but  independent  of  composition.   This  assumption  is  justified 
because  of  the  small  heat  of  mixing  and  the  rather  limited  concentration 
range  of  approximately  from  3.57.  to  77.  encountered  in  this  process. 

A  stage  within  each  effect  consists  of  a  flashing  chamber  and  a  con- 
denser chamber  and  a  demister  which  separates  the  two  chambers.  When  the 
flashing  brine  leaving  one  stage  (the  (n-1)  th  stage)  is  released  into  the 
next  stage  (the  n-th  stage),  water  vapor  flashed  out  of  the  solution.   The 
water  vapor  then  passes  through  the  demister  on  the  way  to  the  condenser 
chamber,  where  it  is  condensed  to  heat  the  nonflashing  brine,  i.e.,  the 
feed  and  recycle  streams. 

3.2   Explanation  of  the  SUBROUTINE  MODEL  Used  in  the  POP-II  Program  ; 

Optimize  the  MEMS  Plant 

Fan,  et.  al.  (7)  have  given  the  theoretical  background  and  the 
derivation  of  the  equations  that  are  given  in  this  section.   The  POP-II 
SUBROUTINE  MODEL  used  in  this  example  is  listed  in  Table  2.   The  equations 
used  in  the  SUBROUTINE  MODEL  are  given  immediately  after  the  listing  of 
SUBROUTINE  MODEL  in  a  one-to-one  correspondence  with  the  FORTRAN  statements 
of  SUBROUTINE  MODEL.   Also,  some  of  the  equation  numbers  are  given  beside 
the  corresponding  SUBROUTINE  MODEL  FORTRAN  statements.   For  example,  the 
first  FORTRAN  statement  of  SUBROUTINE  MODEL  under  the  PERFORMANCE  EQUATIONS 
heading,  is  the  same  as  Equation  (33),  the  first  equation  in  the  following 
explanation  of  SUBROUTINE  MODEL.   It  is  hoped  that  this  method  of  presentation 


2a 


will  give  a  clearer  understanding  of  how  a  SUBROUTINE  MODEL  may  be  programed. 

The  symbols  used  in  SUBROUTINE  MODEL  are  given  in  Table  III-l  of 
Appendix  III.   Because  of  the  simplicity  of  SUBROUTINE  MODEL,  no  logic 
diagram  is  given  for  it. 
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fable   2.    SUBROUTINE   MODEL  .-'or  The   MEMS' Plant 


SUE ROUTINE    MODEL 

COMMON  P(5244) 

DIMFNSICN    X(50),  Y(5C),  CONST(400) 

EQUIVALENCE  (X(1).P(5120),(Y(1),P(5180)) ,( CONST ( 1) ,P(451) ) 

C      SIMULATION  PROGRAM  MEMS  SEAWATER  DISTILLATION  PLANT 

C      ********************  ******  ********  ******************** 

C      CONSTANTS  USED  IN  PERFORMANCE  EQUATIONS 
AL=1000. 
3=17.01723 
AN1=23. 
AN2=23. 
AN3=22. 
CF=0.035 
CP-1. 
XA=230. 
XB=250. 
XC=270. 
XD=290. 
YA=958.8 
YS=945.5 
YC=931-8 
YD=917.5 
U0=510. 

ui=sm. 

■  u  ? = ■>  J  0  . 

U3=510. 

C        ****  I-************************************************* 

C  PERFORMANCE    EQUATIONS 

Y(2)=X(1)+X(  2)+X(3)  (33) 

Y(3)=Y(2)/( l.-CF/X(3) ) 

AL10  =  Y(3)-X(  1)  (35) 

AL2O=AL10-X(2) 

AL3C=AL20-X(3) 

ALl!=Y(3)+X(4> 

AL2I=Y(3)+X( 5)-X(l ) 

AL3I=Y(3)+X(6)-X!1)-X(2) 

CSF10  =  CF*Y!3J/AL10  (41) 

CSF^O=CF*Y(3)/AL20 

CSF1I  =  CSF10*CA:.1C+X(4)  J/AL1I 

CSF2I=CSF2  0*(AL20+X( 5) 1/AL2I 

CSF3I=X( 3)*( AL30+XI6 ) 1/AL3I 

Y{4)=XP.O )-< AL/CP)*AL0G(CSF1C/CSF1I )  (46) 

Y(5)=Y(4)-<AL/CP)*ALGGCCSF20/CSF2I ) 

Y(6)*Yt5)-(AL/CP)*ALCG(X<8 I/CSF3I ) 

A:=:.OIOO+(CSF1I+CSF10)/(2.*0.0300)  (49) 

A2=1.0C75+(CSF2I+CSF2G)/(2.*0.0347> 

A3=0.320l+!CSF3I+X(8) 1/(2. *C. 0315) 

Y(7)=Y(4)-A1  (52) 
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....  2.   (tion't) 


Y(8)=Y(5)-A2 
Y <9)=Y(6)-A3 
2  =  X(-->) 

LAGRANGIAN  POLYNOMIAL 

AA=YA*(Z-XB)*(Z-XC)*(Z-XD)/(  ( XA-XB ) * ( XA-XC ) *  I  X A-XD )  ) 

AB=YR*(Z-XA)*CZ-XC)*(Z-XD)/I ! XB-XA )* ( XB-XC) *< XB-XD) ) 

AC=YC*(Z-XA)*(Z-XB)*(Z-XD)/< 1 XC-XA >* ( XC-XB)*(XC-XD> ) 

AD=YD*(Z-XA>#(Z-XB)*(Z-XC)/( ( XO-XA )* ( XD-XB ) *( XD-XC ) ) 

ALS=AA+AB+AC+AD 

********************************  **  **************  ****** 

YC0)=X(7)*ALS  (55) 

TLF0I=Xt  1C]-Y<  10)/(CP*(Y!3)+XC4>  )  )  (56) 

TLF1I=(CP*( ( AL10+X(5 ) ) *Y (4 )+X t 1 ) *Y ( 7) )-Y( 10 ) ) 
1         /iC?iiV(3]+X(S))l 

TLF2!=(CP*( [AL20+X(6 ) ) *Y ( 5 ) + ( X ! 1 ) +X ( 2 ) ) *Y ( 8 ) ) -Y ( 10 ) ) 
1        /(C?*(Y(3)+X(6) ) ) 

Y(15)=(Y( 10)+Y(3)*CP*X(ll)-(CP*(AL3O*Y(6)  +Y( 2 ) *Y( 9) ) ) ) 
1         /(CP*(Y<6)  -X(ll)  )  )  (59) 

DT0=X(9)-0.5*(X( 101+TLF0I )  (60) 

DT1=X( 10)-TLFCI-Ai-(X(10)-Y(4) )/(2.*ANl ) 

DT2=Y(4)-TLF1I-A2-(Y(4)-Y( 5 ) ) / I 2.*AN2 ) 

DT3=Y(5)-TLr2I-A3-(Y(5)-Y(6)  )/(2.*AN3) 

Y(  11  )=Y(  10)/[DT0*U0)  (64) 

'  Y(12)=X( l)*AL/OTl*Ul) 

Y!13)=X(2)*AL/(DT2*U2) 

Y  tl4)=X(3)#AL/(DT  3*113) 

HP1=X(4)*3*(EXP  t-AL/(0.1104*(X(lC 1+460. )) )  (68) 

1      -  EX?  (-AL/(C.1104*(Y<4)  +46C.)!)) 

HP2=XC5)*B*(EXP  (-AL/(0.1104*(Y(4)  +460.))) 
1      -  EXP  (-AL/(C1104*(Y(5)  +460.))!) 

H?3=X(6)*8*(EXP  (-AL/(0.1104*(  Y(5)  +460.))) 
1      -  EXP  1-AL/(0.1104*(Y(6)  +460.)))) 

****************************************************** 

C0S7  EQUATIONS 

Y(16)=C0NST( 1)*Y(3)/1.E+10  (71) 

Y(17)=CCNST(2)*X(7)/1.E+10 

Y(18)=CO\ST(3)*HP1/1.E+10  (73) 

Y(19)=CCNST(3)*HP2/1.E+10 

Y(20)=COMSTC3)*HP3/l.E+10 

Y ( 2 1 ) =CCN ST { 4 ) * Y ; 1 5 ) / 1 . F+ 1 0 

Y{22)=CCNST(6)*Y(11)/1.E+10 

Y(23)»CCNST(5)«YI  12J/1.E+10  (76) 

Y(24)=CONST(5)*Y< 13  1/1.E+10 

Y(25)=CCNST( 5)*Y(14) /l.E+10 

Y(1)=Y< 16  1+Y1 17i+Y( 18)+Y(19)+Y;20)+Y!21 ) 
1       +Y<22)+Y(23)+YI24)+Y(25)+CONST(7)/l.E-t-10         (77) 

RETURN 
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POP-II   PERFORMANCE  EQUATIONS  FOR  A  MEMS  PLANT 

The  total  production  rate  of  fresh  water,  Cw  is  equal  to  the  sua 
of  the  individual  distilled  water  streams  from  each  effect 


Z\     ■  Wx  +  W2  +  W,  (33) 

A  salt  material  balance"  around  the  entire  plant  gives  the  seaw2ter 
feed  rate 

i  - (34) 

<cf)k 

The  flow  rates  of  the  effluent  flashing  brine  streams  for  the  individual 
effects  are  given,  respectively,  by 

(L)c  =  F  -  WL  (35) 

(L)p  -   (L)c  -  W2  (36) 

(L)K  =   (L)p  -  W3  (37) 

The  flow  rates  of  the  influent  flashing  brine  streams  for  the  individual 
effects  are  given,  respectively,  by 

(L)A  =  F  +  Rj_  (38) 

(L)D  =  F  +  R2  -  Wx  (39) 

(L)G  =  F  +  R3  -  Wj_  -  W2  (40) 


In  the  following  discussion,  the  concentration  of  a  brine  solution  refers 
to  the  salinity  as  defined  by  Badger  and  Associates  (8). 
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Salinity  material  balances  give  the  concentrations  of  the  flashing 

brine  solutions  lea 

ving  effects  one  and  two  respectively: 

C0  F/(L)C                                   (41) 

CQ  F/(L)p                                   (42) 

The  concentrations  for  the  flashing  brine  solutions  entering  effects 

one,  two,  and  three 

,  respectively,  are  given  by  the  following  material 

balance  equations : 

(C  )   = 

A 

(Cf)   (  (L)c  +  RL  )/(L)A                      (43) 

(C-)   (  (L)„  +  R,  )/(L)                       (44) 

t   F      n          J.              D 

(Cf)   (  (L)K  +  R3  )/(L)G                      (45) 

K 

As  the  flashin 

g  brine  proceeds  through  the  effect,  it  is  cooled  due 

to  flashing.   Fan, 

et.  al.  (7)  derived  an  equation  which  gives  the  tempera- 

ture  drop  of  the  fl 

ashing  brine  solution  across  an  effect  as  a  function  of 

the  ratio  of  the  effect's  outlet  to  inlet  concentration.  The  equation  can 

be  rearranged  to  give  the  temperatures  of  the  effluent  flashing  brine 

solutions  for  the  f 

irst,  second,  and  third  effects  respectively: 

(Tf)    - 
C 

(T.)   -   A/C   ln(  (C  )  /(C)   )               (46) 
A        P         C     A 

(T.)   -   A/C   ln(  (C.)  /(CJ   )               (47) 
1  C        P       f  F    £  D 

(Tf)   = 
K 

(T,)   -  A/C     ln(  (C.)  /(C,)      )               (45) 
f  F        p       r  K   £  G 

where  A     is  the  latent  heat  of  evaporation  of  steam  in  the  effect  and  C 

P 

is  the  heat  capacit 

y  of  the  flashing  brine. 

' 
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The  boiling  point  elevation  in  an  effect  is  taken  as  a  linear  function 
of  the  average  flashing  brine  concentration  in  the  effect.   Furthermore,  it 
is  assumed  that  a  one  degree  Fahrenheit  temperature  drop  occurs  across  the 
demister  which  is  placed  in  each  effect  to  prevent  the  flashing  brine 
solution  from  splashing  onto  the  overhead  heat  transfer  area  and  into  the 
distilled  water  collection  troughs. 

Badger  and  Associates  (8)  give  graphs  which  can  be  linearized  about  a 
point  to  give  the  boiling  point  elevation  in  an  effect  as  a  linear  function 
of  the  average  flashing  brine  concentration. 

Equations  (49),  (50),  and  (51)  were  formulated  to  calculate  the  com- 
bined effects  of  the  average  boiling  point  elevation  and  the  demister  temper- 
ature drop  for  effects  one,  two,  and  three  respectively 

(    (Cf)      +  (Cf)    ) 

«1»"    *     l'Ol0O+fos55ff      "\       *  C  <4S> 

(  (c  )    +  (c  )  ) 
<*V    -    whof  ^2         "  (50) 

(    (C   )      +  (C-)    ) 

W)    ,av     =     0.3201  +  " G  r  K  (51) 

0.0315  , 

where  (°0  ,av  is  the  combined  effect  of  the  average  boiling  point  elevation 

and  the  demister  temperature  drop  for  the  n-th  effect. 

The  temperature  of  the  distilled  water  leaving  a  stage  is  given  by  the 

flashing  brine  temperature  leaving  the  stage  minus  the  corresponding  (°c)  ,av 

for  that  stage. 

(T  )   =   (T  )   -  (cC)  ,av  (52) 

C  C       X  C      L 

(T  )  ■>     (T  )  -  (ct)  ay  (53) 

F         F 

(T  )   -  (T  )  -  (<*),,av  (54) 

K       r  K      J 
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The  latent  heat  of  vaporization  of  the  saturated  brine  heater  steam 
is  written  as  a  function  of  the  steam  temperature  by  using  a  Lagrangian 
Polynomial.   However,  in  this  study  the  brine  heater  steam  temperature  is 
allowed  to  vary  only  in  one  case  study  and  then  by  only  a  small  amount, 
because  increasing  the  temperature  will  bring  about  a  significant  increase 
in  the  pressure  within  the  brine  heater.  Higher  pressures  require  that  the 
brine  heater  be  built  out  of  sturdier,  more  expensive  material.  The 
effect  of  brine  heater  pressure  upon  brine  heater  cost  has  not  been  taken 
into  account  except  at  the  brine  heater  temperature  used  for  this  study. 

The  rate  of  heat  addition  to  the  brine  heater  by  the  brine  heater 

steam  is  given  as 

q   =  m*  (55) 

^s      s  s 

where  m  is  the  steam  flow  rate  and  X   is  the  latent  heat  of  vaporization 

of  the  brine  heater  steam. 

A  heat  balance  around  the  brine  heater,  including  all  effects  up  to 

the  (n  +  l)th  effect,  will  give  the  temperature  of  the  seawater  and  recycle 

brine  entering  the  n-th  effect.   Equations  (56),  (57),  and  (58)  give  the 

temperatures  of  the  combined  non-flashing  brine  streams  entering  the  brine 

heater,   effect  one,  and  effect  two,  respectively.   The  seawater  coming 

into  the  plant  is  assumed  to  be  at  an  ambient  seawater  temperature  of  85°  F. 

(T. )   =   (T.) 


J 'A       fA    Cp(F+V  (56) 

C   (((L  )   +R  )  (T  )   +  W  (T  )  )  -  q 
(T.)    =   P      t  C     Z     r  C     X   C  C  (57) 

2  c  CP  <F  +  V 

C   (((L  )  +  R  )  (T  )  +  (W  +W,)(T  )  -a 
(T.)   -   P     t   F    J    r  F     l  -   I   F  '      (53) 

F  C   (F  +  R  ) 

P       -1 
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A  heat  balance  around  the  entire  plant  will  give  the  required  cooling 
water  flow  rate 


<!„  +  FC„  <T 


)   -  C   f (L)  (T,)  +  (  W  )  (T  )  ] 
J  K    p  I         K  t  n     c  RJ 


R4  =         I        J  K    '■"I.    2—1 -  :;      <59) 

(TJ    -   (Ts)  I 


f  k  j  k; 

The  temperature  driving  force  available  for  heat  transfer  in  the 
brine  heater  is  taken  as  the  arithmetic  mean  of  the  approach  temperature 
differences. 

(At)Q  =  T=  -(1/2)((T.)   +  (T  )   )  (60) 

A       •'A 

Fan,  et.al.  (7)  derive  equations  for  calculating  the  effective  heat 

transfer  driving  force  for  each  effect.  The  equations  take  into  account 

the  boiling  point  elevation  of  the  flashing  brine,  the  demister  temperature 

drop,  and  the  fact  that  the  effect  has  a  finite  number  of  stages  rather  than 

an  infinite  number.  The  effective  heat  transfer  driving  forces  for  effects 

one,  two,  and  three  are  given  by  Equations  (61),  (62),  and  (63)  respectively 

(Tf)   -  (Tf) 

(At)   =   (T.)   -  (T  )   -  (oC)   av  -(1/2)— i— — L    (61) 

A     J  A  '1 

(Tf)   -  (Tf) 
(At),  -   (Tf)   -  (T.)   -(<*),, av  -(1/2) 2 L.    (62) 

»     JB  N2 

<Tf)   -  (Tf) 
(At)   -   (T  )   -  (T  )   -  (cC)   av  -(l/2)__£ _L    (63) 

where,  N  is  the  number  of  stages  in  effect  n. 

The  heat  transfer  areas  required  for  the  brine  heater  and  each  effect 
are  given  as  follows 

A0  -  qs  /  (At)Q  UQ  (64) 

Ax  =  Vx     /    (flt)L  UL  (65) 
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A2  =  W2  /  (At)2  U2  (66) 

A3  »  W3  /  (At)3  U3  (67) 

where  U  is  the  overall  heat  transfer  coefficient  for  the  brine  heater  at 
n  th  effect.   Equations  (65),  (66),  and  (67)  give  a  somewhat  approximate 
evaluation  of  the  heat  transfer  areas  of  the  heat  rejecting  sections  of 
the  respective  effects.  This  approximation  is  thought  to  be  justified 
since  the  sizes  and  thus  the  costs  of  the  heat  rejection  areas  are  a  small 
part  of  the  total  heat  transfer  area  of  each  effect. 

The  power  requirement  of  a  recycle  brine  pump  is  taken  to  be  propor- 
tional to  the  product  of  the  pressure  head  which  is  developed  by  the  pump 
and  the  recycle  brine  mass  flow  rate.  The  pressure  head  developed  by  the 
pump  must  be  enough  to  overcome  the  frictional  loss  in  the  recirculation 
line  and  the  pressure  increase  caused  by  different  vapor  pressures  at  the 
inlet  and  outlet  of  the  stage.  In  order  to  determine  the  power  requirements, 
the  following  assumptions  are  made: 

1.  The  vapor  pressure  of  the  flashing  brine  solution  is  given 
by  an  integrated  Clausius-Clapeyron  equation. 

2.  The  friction  loss  is  proportional  to  the  pressure  drop 
between  the  outlet  and  inlet  of  the  effect. 

3.  The  density  of  the  flashing  brine  solutions  are  constant 
throughout  the  plant. 

4.  The  horsepower  calculated  has  been  corrected  for  the 
mechanical  inefficiency  of  the  pump. 

The  power  requirements  for  effects  one,  two,  and  three,  respectively,  are 
given  oy 
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(HP), 


"1   "1 


Exp 


0.1104  (  (Tf)  +  460.) 

A 


(HP)2  =  Bx  R2 


(HP)3  -  B]_  R3 


-Exp 


Exp 


-Exp 


Exp 


-Exp 


-\ 


0.1104  (  (T  )  +  460.) 
C 


0.1104  (  (T  )  +  460.) 
C 


0.1104  (  (Tj  +  460.) 


0.1104  (  (T.)  +  460.) 

F 


0.1104  (  (T  )  +  460.) 
K 


(68) 


(69) 


(70) 


where  (HP)  is  the  power  requirement  for  effect  n  and  B.  is  a  constant 
which  takes  into  account  the  friction  loss,  the  recirculating  brine  density, 
the  mechanical  efficiency  of  the  pump,  and  the  integration  constant  of  the 
integrated  Clausius-Clapeyron  equation. 


POP-II  ECONOMIC  EQUATIONS  FOR  A  MEMS  PLANT 

Both  initial  equipment  costs  and  operating  costs  for  the  life  of  the 
equipment  are  taken  into  account  in  the  objective  function.  The  following 
costs  are  considered  to  be  significant  in  the  selection  of  equipment  size: 
I.   Initial  equipment  costs 

a.  Heat  transfer  area  costs 

1.  Brine  heater 

2.  Each  Effect 

b.  Pumps  -  3rine  recirculation 

c.  Outer  shell  of  each  effect 
II.   Operating  costs 

a.  Feed  brine  pretreatment  and  pumping 

b.  Cooling  water  pumping 

c.  Recirculation  brine  pumping 

d.  Brine  heater  steam 

In  the  following  cost  equations,  the  cost  coefficients  have  been 
divided  by  the  total  production  rate  of  fresh  water  to  obtain  the  total 
cost  on  the  basis  of  a  unit  of  production. 

The  initial  cost  of  the  brine  feed  pump,  the  brine  pretreatment  cost,' 
and  the  brine  feed  pump  operating  cost  are  all  considered  to  be  proportional 
to  the  brine  feed  rate.  Therefore, 

Ei  =  V  CD 

where  C_,  is  the  unit  feed  water  cost. 

w 

The  brine  heater  steam  cost  is  proportional  to  the  steam  consumption 

rate 

E2  -  Csms  <72> 

where  C  is  the  unit  steam  cost. 
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The 

energy  cost  associated  with  the  operation  of  the 

recycle  purr.p  in 

the  n-th 

effect  and  the  initial  cost  of  each  pump  can  be 

calculated  by  the 

following  equation: 

(ce  +  YCj)  (hp)b  =  cHp  (IIP)  n 

(73) 

where  C 
e 

is  the  unit  power  cost,  r  is  the  capitalization 

charge,  C  is  the 

capital  cost  per  horsepower,  and  C   =   (C  +  ^  C  ). 

The 

cooling  water  cost  is  directly  proportional  to  the  quantity  used 

E,  =  Cr  R. 

3      u  4 

(74) 

where  C„ 

is  the  unit  cooling  water  cost. 

The 

brine  heater  capital  cost  per  hour,  E«,  is  given 

by 

E2   -  VCBAQ   =   CktAq 

(75) 

where  C_ 

is  the  capital  cost  per  unit  heat  transfer  area, 

andCHT=yV 

The 

heat  transfer  area  cost  for  the  n-th  effect  can  be  calculated  as 

follows: 

<EV  =YVn  =  CeA 

(76) 

where  A 

n 

is  the  heat  transfer  area  of  the  n-th  effect,  C„ 

is  the  capital 

cost  per 

unit  heat  transfer  area,  and  C__  =  j   C  .  There  should  be  a 

different  cost  coefficient  for  the  brine  heater  area  cost 

than  for  the 

effect  areas. 

The 

total  outer  shell  cost  for  the  plant  based  on  a  unit  of  production 

is  taken 

into  account  by  adding  a  constant,  E, ,  the  sum  o] 

the  previous 

costs.  The  total  water  cost  per  unit  of  production  is  given  by  the  sum 

of  the  cost  terms  previously  stated. 
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The  objective  cost  to  be  minimized  is: 

C„F  +  C.m_  +   g  CHp(KP)n  +  CcR4  +  C^  + 


n+1 


N=3 

C  c  A  +  E„  (77) 

n=l  EF  n    u 


3.3  Results,  Discussion  and  Conclusions  Concerning  the  Usefulness  of 
POP-II  for  Optimizing  the  MEMS  Plant 

The  results  of  three  representative  runs  made  with  the  POP-II 
optimization  program  are  summarized  in  Table  3.  Runs  1,  2,  and  3  were 
started  at  different  feasible  points.   The  initial  and  final  feasible 
X(J),  values  and  the  total  cost  obtained  for  each  run  are  given  in  the 
table.   Three  independent  variables,  the  brine  heater  steam  temperature, 
X(9),  the  temperature  of  the  flashing  brine  entering  effect  one,  X(10), 
and  the  temperature  of  the  fresh  seawater  entering  the  plant,  X(ll),  were 
held  constant  in  this  work  so  tha..  the'  results  could  be  compared  with  the 
work  by  Fan,  et.  al.  (7).   The  complete  computer  output  for  Run  1,  includ- 
ing the  loop-to-loop  output  is  given  in  Appendix  II.   The  loop-to-loop 
outputs  for  Runs  2  and  3  are  also  given  in  Appendix  II. 

The  starting  values  of  X(J)  used  in  Run  1  were  essentially  the  same 
as  the  overall  optimum  for  the  MEMS  process  obtained  by  Fan,  et.  al.  (7) 
using  the  discrete  version  of  the  maximum  principle.  The  total  cost  for 
this  set  of  initial  input  values  was  ?0.2867  per  thousand  gallons  of  potable 
water,  which  is  again  very  close  to  the  optimal  value  obtained  by  Fan, 
et.  al.  As  shown  in  Table  3,  an  appreciable  change  occurred  between  the 
initial  and  final  values  of  only  two  independent  variables,  X(7)=  m  and 
X(8)=(Cf)k.   The  brine  heater  steam  consumption  rate,  m  ,  decreased  by 
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1.9  percent.   The  concentration  of  the  discharge  brine  solution,  (C_)  , 
increased  by  14.1  percent.   The  final  total  cost  was  $0.2866  per  thousand 
gallons  of  potable  water  produced  which  is  again  very  close  but  slightly 
higher  than  that  obtained  by  Fan,  et.  al.   Since  Run  1  was  the  lowest  cost 
case  obtained  with  POP,  the  final  values  of  the  dependent  variables  are 
given  in  Table  4. 

The  total  cost  essentially  stayed  constant  throughout  Run  1.  No 
infeasible  intermediate  loops  were  encountered  during  the  run  as  can  be 
seen  by  noting  that  the  production  equality  constraint  was  always  satisfied 
and  the  independent  and  dependent  variables  remained  positive  throughout 
the  run.   The  adaptive  move  limit  option  was  used  in  this  run.   The  run 
terminated  at  the  end  of  twenty-five  loops  because  the  calculated  move 
limits  became  too  small. 

Run  2  was  started  with  all  X(J)  values  except  X(3)  less  than  the 
corresponding  starting  X(J)  values  of  Run  1.   The  initial  total  cost  was 
$0.2926  per  thousand  gallons  of  potable  water  produced.  As  was  the  case 
in  Run  1,  only  the  brine  heater  steam  consumption  rate  and  the  concentration 
of  the  discharge  brine  solution  changed  appreciably  during  the  run.   The 
brine  heater  steam  consumption  rate  increased  by  10.8  percent.   The  con- 
centration of  the  discharge  brine  solution  increased  by  14.6  percent.   The 
final  total  cost,  $0.2890  par  thousand  gallons  of  potable  water  produced, 
was  only  0.84  percent  larger  than  the  optimum  cost  found  in  Run  1. 

During  the  course  of  Run  2,  several  loops  ended  at  infeasible  points. 
For  these  loops,  the  process  required  a  negative  amount  of  cooling  water  in 
effect  three,  R4#   The  next  loop  after  an  infeasible  loop  was  always  feasible 
and  the  last  ten  loops  were  all  feasible.  Also,  small  oscillations  in  the 
total  cost  in  going  from  one  loop  to  the  next  were  noted  at  the  beginning 
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of  the  run.   This  run  also  terminated  at  the  end  of  twenty-five  loops 
because  the  calculated  move  limits  became  too  small. 

Run  3  was  started  with  all  X(J)  values  except  X(3)  greater  than  the 
corresponding  X(J)  values  of  Run  1.   The  initial  cost  was  $0.3340  per  thous- 
and gallons  of  potable  water  produced.  In  this  run  all  of  the  independent 
variables  changed  appreciably  during  the  run.   The  significant  changes 
in  the  values  of  the  independent  variables  can  perhaps  be  explained  by  the 
fact  that  the  error  limits  on  some  of  the  dependent  variables  were  larger 
in  this  run  than  in  Runs  1  and  2.   Furthermore,  the  initial  move  limits  of 
the  independent  variables  were  larger  in  this  run  than  in  Runs  1  and  2. 
Loop  number  six  was  the  best  loop  obtained  in  this  run.  The  total  cost 
at  the  end  of  loop  six  was  $0,306  per  thousand  gallons  of  potable  water 
produced.   Even  though  this  run  ended  infeasible,  the  total  cost  at  the 
end  of  loop  six  was  only  6.S  percent  higher  than  the  optimum  cost  obtained 
in  Run  1. 

At  the  end  of  loop  one  of  Run  3,  the  total  cost  was  nearly  three  times 
the  starting  total  cost.   Thereafter,  the  total  cost  progressively  decreased 
until  the  end  of  loop  six.   The  infeasible  loops  of  this  run  had  negative 
cooling  water  requirements  for  effect  three.   This  run  terminated  at  the 
end  of  eleven  loops  because  the  calculated  move  limits  became  too  small. 

Preliminary  simulation  of  the  operation  of  the  MEMS  plant  was  carried 
out  using  computer  program  MODEL.   The  main  portion  of  MODEL  consists  of 
the  FORTRAN  statements  used  in  SUBROUTINE  MODEL.  MODEL'S  main  function 
was  to  debug  the  FORTRAN  statements  used  in  SUBROUTINE  MODEL.   Since  remote 
computer  facilities  were  used  to  run  the  POP  program,  it  was  decided  to 
debug  SUBROUTINE  MODEL'S  FORTRAN  statements  on  available  computer  facilities 
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in  order  to  speed  up  the  debugging  process.  Normally,  the  POP  program 
may  be  used  to  do  the  debugging.  The  starting  points  of  Runs  2  and  3 
were  also  determined  by  using  MODEL.   MODEL  is  explained  in  Appendix  III. 

The  loop-to-loop  output  is  a  very  useful  aid  for  studying  a  particular 
SUBROUTINE  MODEL.  Three  general  observations  about  the  MEMS  SUBROUTINE 
MODEL  were  made  by  inspecting  the  loop-to-loop  outputs  for  the  three  runs. 
First,  there  appeared  to  be  a  strong  interaction  between  the  brine  heater 
steam  consumption  rate,  m  ,  and  the  discharge  brine  concentration,  (C-),  . 
This  was  readily  seen  by  observing  the  loop-to-loop  output  of  Run  2. 
Secondly,  small  charges  in  m  and  (C_).  in  m  and  (C„).  caused  rather 
large  changes  in  the  cooling  water  requirement,  R, .  Most  of  the  runs 
which  became  infeasible  did  so  because  R,  became  negative.  Many  tines 
R,  became  negative  after  only  a  small  change  in  either  m  or  (Cf),  . 

The  DYDX  matrix  is  also  very  useful  for  studying  a  particular  SUB- 
ROUTINE MODEL.   For  the  MEMS  SUBROUTINE  MODEL  most  of  the  a.,  elements  of 
the  m  and  (C.),  columns  of  the  DYDX  matrices  of  all  three  runs  were  large 
in  comparison  to  the  a,,  elements  of  the  other  independent  variables. 
Consequently,  as  was  observed  in  the  loop-to-loop  outputs,  small  changes  in 
these  independent  variables  caused  relatively  large  changes  in  the  other 

dependent  variables  in  addition  to  the  total  cost.  Most  of  the  a.,  elercnts 

ij 

of  the  cooling  water,  R, ,  row  of  the  DYDX  matrices  were  larger  than  the 
elements  of  the  other  rows.  This  explains  why  small  changes  in  most  of  the 
independent  variables  had  a  large  effect  on  the  cooling  water  requirements. 

Several  conclusions  regarding  the  usefulness  of  POP  have  been  formed 
in  working  with  the  MEMS  plant  problem.   It  is  not  a  routine  matter  to 
apply  POP  to  optimize  complex  processes.   Experienced  programmers  may  even 
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be  required  to  use  POP  on  simple  models.  It  appears  that  some  experimenta- 
tion with  the  modes  of  operation  and  the  input  data  may  be  required  on 
simple  models. 

In  the  next  section,  future  work  is  proposed  in  which  changes  in  the 
mode  of  operation  and  the  input  data  cards  can  be  made  now  that  more  infor- 
mation is  available  on  the  SUBROUTINE  MODEL  used  in  this  work. 

3.4  Suggestions  for  Further  Work 

There  are  several  possible  changes  that  can  be  made  in  the  POP  input 
data  forms  and  in  the  manner  of  operating  POP  which  may  help  POP  to  con- 
verge to  the  optimum  when  it  is  started  at  points  which  are  far  away  from 
the  optimum. 

A  minimum  limit  can  be  put  on  the  dependent  variable  R,  which  hopefully 
will  prevent  it  from  becoming  negative.  The  optimizer  may  be  run  without 
using  the  adaptive  move  limit  calculation  feature.  This  mode  of  operation 
should  enable  the  optimizer  to  get  closer  to  the  minimum  cost  during  a  run. 
Since  the  model  is  sensitive  to  small  changes  in  the  independent  variable 
(C„),  ,  this  variable  can  be  held  at  a  different  constant  value  for  a  series 
of  runs  as  was  done  by  Fan  et.  al.  (7).   Furthermore,  a  reformulation  of 
SUBROUTINE  MODEL  using  different  independent  variables  may  be  advisable. 

Oscillations  in  the  total  cost  may  possibly  be  reduced  by  making  the 
following  adjustments:  Move  limits  on  the  less  sensitive  independent  vari- 
ables may  be  increased.   It  may  be  necessary  to  individually  adjust  each 
move  limit  during  a  series  of  runs.   In  combination  with  the  above  steps, 
changes  in  the  values  of  the  constants  used  in  the  optimizer's  convergence 
logic  may  help. 

A  general  procedure  for  conducting  POP  runs  is  outlined  next.   Analyze 
the  loop-to-loop  output  and  the  DYDX  matrix.  Make  necessary  adjustments  such 
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as  those  previously  described.  Restart  the  next  run  from  the  best  point 
obtained  in  the  previous  run.  Many  runs  which  do  not  end  at  a  feasible 
point  contain  useful  information  and  often  one  of  their  loops  is  an  ex- 
cellent starting  point  for  the  next  run.   See  loop  six  of  Run  3.   Sample 
input  data  forms,  which  contain  the  above  changes  are  given  in  Appendix  I. 

A  more  elaborate  process  simulation  program  or  SUBROUTINE  MODEL  can 
be  investigated  in  the  future.   The  number  of  stages  per  effect  can  easily 
be  made  an  independent  variable  as  is  suggested  by  Fan,  et.  al.  (7). 
Information  giving  equipment  cost  as  a  function  of  equipment  size  and  steam 
cost  as  a  function  of  the  steam's  physical  properties  (9)  is  becoming  more 
readily  available  as  time  progresses. 

Certain  advantages  may  be  realized  by  combining  POP  with  other  optimi- 
zation techniques  such  as  the  discrete  maximum  principle. 
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5.0  NOMENCLATURE 

ai         A  constant  in  several  equations. 
eij         As  used  in  equations  (15)  and  (17),  a  constant 
which  denotes  the  coefficient  of  the  i  th  row 
and  the  j  th  column  of  a  system  of  linear 
equations.  As  used  in  equation  (26)  it  repre- 
sents the  partial  derivative  of  the  i  th 
equality  constraint  with  respect  to  the   j  th 
independent  variable. 
^0  Keat  transfer  area  of  the  brine  heater,  sq.  ft. 

Al  Heat  transfer  area  of  the  first  effect,  sq.  ft. 

A2  --eat  transfer  area  of  the  second  effect,  sq.  ft. 

A3         H*at  transfer  area  of  the  third  effect,  sq.  ft. 
bi         A  constant  which  is  first  introduced  in  the 

equality  constraint,  equation  (lb). 
Bl         A  constant  used  in  the  equations  which  calculate 
the  horsepower  requirements  of  the  recycle  pumps, 
H.?.-hr./lb. 
ce         Unit  power  cost,  $/H.P. 
(Cf)A       Concentration  of  the  flashing  brine  stream  flowing 

into  the  first  effect,  wt,  %. 
(-f )c       Concentration  of  the  flashing  brine  stream  flowing 
out  of  the  first  effect,  wt,  %. 


<Cf>D 

<cf)F 

<Cf>G 

<cf>K 

Cij 

co 

S 

cs 

CB 

CC               ' 

C^,-               : 

GH               , 

- 

■      CHT             I 

CJ            c 

°w           l 

- 

Goncentratior.  of 

the  flashing 

brine  stream 

flowing  into  the 

second  effect,  wt.  %, 

Concentration  of 

the  flashing 

brine  stream 

flowing  from  the 

second  effect,  wt.  %, 

Concentration  of 

the  flashing 

brine  stream 

flowing  into  the 

third  effect, 

wt.  %. 

Concentration  of 

the  f lashing 

brine  stream 

flowing  from  the 

third  effect, 

wt.  %. 

The  coefficient  of  the  quadratic  tern  of 

equation  (16) 

Concentration  of 

seawater  fed 

to  the  MEMS 

Plant,  wt.  %. 

Heat  capacity  of 

flashing  brine  solution,  3tu/ 

lb.  °P. 

Unit  steam  cost, 

?/lb. 

Unit  cost  of  brir 

e  heater  heat 

transfer  area, 

S/sq.  ft. 

Unit  cooling  water  cost,  $/lb. 

Equal  to  Y  CH 

Unit  cost  of  heat 

transfer  are 

a  in  an  effect, 

?/sq.  ft. 

Equal  to  (Ce+fCj 

) 

Equal  to  YC-, 

Capital  cost  per 

horsepower,  $ 

/:-:.?. 

Unit  feed  water  cost,  $/lb. 
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CONST(I)         The  i  th  constant  in  the  ^rray  of  constants 

which  can  be  used  with  POP-II 
J*-  -  rrie  name  of  the  matrix  printed  out  at  the  end  of 

a  POP-II  run  which  can  be  used  to  make  an 
incremental  sensativity  analysis  near  the  point 
at  which  the  run  ends. 
ei  Tke  value  of  the  i  th  equality  constraint  when  It 

evaluated  at  a  feasible  point 
El  Total  cost  associated  with  the  seawater  fed  to 

the  MEMS  plant,  5 /1000  gal.  fresh  water 
E2  Total  cost  of  brine  heater  steam,  S/1000  gal. 

fresh  water 
E3  Total  cost  of  cooling  water,  S/1000  gal.  fresh 

water 
u3  ;  Total  cost  of  heat  transfer  area  for  the  a  th  effect 

S/1000  gal.  fresh  water 

Total  cost  of  the  outer  shell  of  the  MEMS  plant 
based  on  1000  gal.  fresh  water  produced, 
3/1000  gal.  fresh  water 
^x'  The  objective  function 

-Jeed  rate  to  the  MEMS  plant,  lbs/hr. 
S±(.x)  The  i  th  equality  or  inequality  constraint 

h^  A  function  used  in  the  discussion  about  the 

quadratic  form 
\(x'  The  *•  th  equality  constraint  (used  in  the 

explanation  of  POP-II) 


-6 


- 


55 

(HP), 

Horsepower  required  for  the  recycle  pump  of  the 

first  effect,  H.?. 

(HP)2 

Horsepower  required  for  the  recycle  pump  of  the 

second  effect,  H.P. 

(H?)3 

Horsepower  required  for  the  recycle  pump  of  the 

third  effect,  H.P. 

L 

The  Lagrangian  function  used  in  the  development 

of  the  Kuhn- Tucker  conditions. 

(L) 

A 

Flow  rate  of  lashing  brine  stream  into  the 

first  effect,  Ibs./hr. 

U)c 

Flow  rate  of  lashing  brine  stream  from  the 

first  effect,  lbs./hr. 

(L)D 

Flow  rate  of  lashing  brine  stream  into  the 

second  effect,  lbs./hr. 

(L)? 

Flow  rate  of  flashing  brine  stream  from  the  second 

effect,  lbs./hr. 

a>G 

Flow  rate  of  flashing  brine  stream  into  the 

third  effect,  lbs./hr. 

<L)K 

Flow  rate  of  flashing  brine  stream  from  the 

third  effect,  lbs/hr. 

m 

s 

Brine  heater  steam  flow  rate,  lbs/hr 

^s 

Heat  transfer  rate  ir.  the  brine  heater,  BTU/hr. 

-_ 

Recycle  brine  flow  rate  in  the  first  effect, 

lbs. /hr . 
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R2 

Recycle  brine  flow  rate  in  the  second  effect, 

lbs./hr. 

Recycle  brine  flow  rate  in  the  third  effect, 

lbs. /hr . 

H 

Cooling  water  flow  rate  in  the  third  effect, 

lbs./hr. 

s 

Value  of  the  objective  function  of  an 

optimization  problem 

(Xc)r 

Temperature  of  distillate  water  from  effect  one, 

°P 

Cfc), 

Temperature  of  distillate  water  from  effect  two, 
Temperature  of  distillate  water  from  effect 

(Tc) 
K 

three,  °F 

*  A 

Temperature  of  flashing  brine  stream  into 

effect  one,  °F 

C«£>q 

Temperature  of  flashing  brine  stream  from  effect 

one,  °F 

(Tf)F 

Temperature  of  flashing  brine  stream  from 

TOO,  °? 

<*f>K 

Temperature  of  flashing  brine  stream  from 

effect  three,  °'S 

<Va 

Temperature  of  combined  recycle  and  feed  streams 

into  the  brine  heater,  °F 

(Tj>c 

Temperature  of  combined  recycle  and  feed  streams 

into  effect  one,  °?  ■ 

„ 

(tj), 

Temperature  of  combined  recycle  and  feed  streams 

into  effect  two,   P 

<VK 

Temperature  of  combined  recycle  and  feed  streams 

into  effect  three,   F 

-s 

Temperature  of  brine  heater  steam,  °F 

(At)0 

Temperature  differences  used  for  heat  transfer 

in  brine  heater,  °? 

CAt)x 

Temperature  difference  used  for  heat  transfer 

in  effect  one,  °? 

(At)2 

Temperature  difference  used  for  heat  transfer  in 

effect  two,  °? 

(At)3 

Temperature  difference  used  for  heat  transfer  in 

effect  three,  °f 

ui 

Slack  variable  used  in  development  of  Kuhn- 

Tucker  Conditions 

U0 

Overall  heat  transfer  coefficient  for  the  brine 

heater,  Btu/hr-f t2_°j? 

°1 

Overall  heat  transfer  coefficient  for  effect  one, 

Btu/hr-ft  -°F 

U2 

Overall  heat  transfer  coefficient  for  effect  two, 

Btu/hr-ft2-°y 

s 

Overall  heat  transfer  coefficient  for  effect 

three,  Btu/hr-f t2-°x'' 

Ul 

Equal  to  ^ 

*1 

Kate  of  distillate  water  production  in  effect  one, 

Ibs./hr. 

£*n 


Rate  of  distillate  water  production  in  effect 

two,  los./.'.r. 

...  te  of  cistillate  water  production  in  effect 

three,  lbs./hr. 

Total  rate  of  distillate  water  production  in 

the  MEMS  plant,  lbs./hr. 
x  An  abbreviated  notation  which  stands  for  the 

point  (x,,  X2 ,  ••••,  xs)  in  s-space. 

xQ  A  feasible  point  (xLq,  x?q,  ,  xgQ) 

X(I)  The  i  th  independent  variable  (used  in  POP 

programming) 
y«  The  i  th  dependent  variable 

yiCx.-r  §x)       The  i  th  dependent  variable  approximated  as  a 

linear  function  of  £x 
Y(l)  The  i  th  dependent  variable  (used  in  POP 

programming ) 
Y(l)  The  objective  function  (used  in  POP  programming) 

Greek  Letters 


(  d.  }  The  average  boiling  point  elevation  in  the  n  th 

r.,av 


effect 

Latent  heat  of  vaporization  of  flashing  brine 
stream  in  an  effect,  otu/lb. 
Lagrange  multiplier 
X  s  Latent  heat  of  steam  used  in  brine  heater, 

Btu/lb. 
Capitalization  charge 


X 


r 
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Ar  INDIX  11 

POI  -XI  3  ...  __  OUTPUT 

Ihe  computer  o t  for  Jun  1  Is       in  this  Ape 

....  ..  _.  ..:.  optimum  cas  -  for  - —  ork.   Also,  the  loop-to- 

loop  computer  outputs  for  Suns  2  and  3  are  given  in  tnis 

..dix. 
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RUN  Ml   FR   1       PROBLEM  NUMBER   !       I  OOP  NUMBER   0 
INPUT  r  y*  I  /»  FOR  PLANT   IPTIKIZATION  PROGRAM,  PAGE    I 


12/12/66 


*****   KU 

-  o  u;  c  -  o  -  o  -  o  -  o   3 

-0    I   -0    1   -0   -i)   -0 


12   12   66   -0    I   -0   -0   -0    9    1    I    I   -0 

-0    5  _-0   -0   -o   -o   -0   -0   -0.  -0    I   -0   -0 


«.**   CS1  -  CPTIMIZER  CCNSTANJS  *****_**J* 

O.i-  994E-04  0.099999994E-04  0.5G0000000E    03 

-  0  .  0  -  1999 9 9 9 961      C  L      -  0  . 


0.9999999R5E-0  2 

-0. 


0-499999993E-0 I 

0.50Q000000E  00 


0.499999993C-02 

o.  n ooooooo or  no 


******  c. c  h s  r  -  mxcbl  0:1. i  r  ic  ii.-k  rs**** 

0.17  64  '  I'  II  0  5    0.25 00 000 00E  0  7    0.575899996E  08 

0  .20  /'.''•  )9  n     09 


0.598399997E    04  0.239699997E    06 


0.375999995F     06 
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IvLN    MNEI  R       l 

INPU r    I.  Al  /<    FOR 


:     I     '  '    ■  Bl  k 
P L ^ i\  I    OPTIMZATIGN 


3  LOOP    NUMBER 

PROGRAM,    PAGI  2 


12/12/66 


:  ■  ■  ■  [NCI  'l  I        LEJ  •■  *  * 


V  AR  I A  E  L I  I  •  ■  II 
.  I 
W2 
V-  ! 
R  I 
R2 

MS 

CSI   1C 
TS 

rsi  ii 

i  ii  n 


>****       CI  It  iv  Dl  N  I     VAP  I  ^  0  L  E  S 


V  A  R  I A  E 1 1     P  \  i  I 
r  C  T  A I     ii  S I 

P        I   II    ll(f 


I  I  I 

TM 
l  SI 
rsi 


D 
ii 

i 
'i 


rcio 
rc3i 

Q 

ari  a 

ARI  A 

■     i  A 
AKI  A 

(ML    Vn 

F    CCS1 

.    ,  i  «  i 

i  i  ■   i 

I 
i  i  •  I 

I 


Rl 
R2 

CW 
AO 


IN 
2.E.6 

8  3  4 

204 

1  5  '< 
10  3 
20  1 
152 
101 
49  1 
J  3  / 
0.992 
0.IW 
0.  1  29 
0,.e76 


i;ur 

7  1L-0  0 
OOE     04 


04 


06£ 

.'■  i 
42J 

: 

ii  n: 

15E 

•  n 
i  /r 

4.CE 

><  i  - 


05 
03 
03 
03 
03 
03 
03 
06 
02 
0  3 
04 
0  4 
04 
01 
0  0 


02 


■02 
02 

-02 
731  -02 
9  61 


MIN,    LIMIT 
•0. 

0  .  ii  3  4  0  0  E 
•0. 
-0. 
■0. 

0 .  8  5  0  0  1  E 
■0. 
■0. 
■0. 
■0. 
■0. 
■0. 
■0. 
•0. 
•0. 
■0. 
•0. 
■0. 
■o. 
•0. 
•0. 
•0. 


MAX. 

0.900 
0.900 
0.  100 
0.900 
0 .  9  0  0 
0.900 
0.900 
0.100 
0.279 
0.255 
0,85  0 


L  i  m  i  r 

OOE    05 


00C 

OOE 

001 

OOE 

OOE 

OOE 

OOE 

4  IF 

OLE 

10E 


05 
0  5 
05 
05 
05 
0  5 
01 
0  3 
03 
02 


MAX.     LIMIT 
-0. 

0.83400E    04 
■0. 

0.  27940 E     0  3 
■0. 
•0. 

0. 

0. 

0. 
•0. 

0. 
•0. 
■0. 
•0. 
■0. 

0. 

0. 

0. 

0. 

0. 

0. 

0. 


0.100 
0.  100 
0.  100 
0.300 
0.300 
0.300 
0  .500 
0.50  0 
0.500 
0.500 
0.500 


DELTA 

OOE    02 


OOE 
OOE 

OOE 

OOE 

OOE 

OOE 

0  0  E-  0  3 

00E-03 

00E-01 

OOE-03 


02 
02 
02 
02 
02 
01 


OELTA 

0 .10000E  01 
0.  10000E  01 
0. 10000E  01 
0.1  00 OOE  0  1 
0.  10000E  0  l 
0.  10000E  01 
0  .  1  0  0 0  0  E  0  1 
0  .  1  0  0  0  0  E  0  1 
0  .  1  0  0  0  0  F  0  1 
0 .  1 0  0  0  0  E  01 
0  .  I  0  0 0 0  F  01 
0 .  1  0  0  0 0  F  01 
0  ,  1  0  0  0  0  E  01 
0.10 0 0 0 E  0 1 
0  .  1  0  0  0  0  F  0  1 
0.10000E  0  1 
0.1 00 OOE  01 
0.1 00 OOE  01 
0.1  0  0  OOF  0  1 
0.10 0 0 0 F  0 1 
0.100  OOF  0  1 
0.10000E    0  1 


IVE 

0.  1 

0.1000 
0.  1000 
0.3000 
0.3000 
0 ..  3  0  0  0 
10  0  0 
5  )  0  0 
5  3  0  0 
5)00 
50  0  0 


L  I  M  I  T 

OF    0  3 


OE 

OF 
OF 
OE 
OF 
OF 


03 
0  3 
03 
03 
0  3 
0  2 


OE-02 
OE-02 
OF-fl? 

OF-02 


ERROR  LIMIT 

1 .  0  0  0  0  0  E-  0  3 
1.00  0  0  0E-03 
1 .  0  0  0  0  0  C-  0  3 
0.50 0  0  OF  0  0 
0.50  0 OOE  0  0 
0.50000E  00 
0.50 0  0  or  0  0 
0.50  0 OOE  0  0 
0.50 0  0 0  E  00 
0.10 0 0  0  E  03 
0.1)0  0  0  '"-00 
0.10  000E  01 

o.noou   oi 

0 .  1 0  0  0  0  E  0  1 
0.100 0  OE    01 

1.  00  00  OE- 0  3 
1.00  0001  03 
1. 00  M 
1.0)  1.3 
I  .  0  0  0  0  0  I  -  0  4 
I. 00 00 Of  -0  ! 
L.OOOOOE-04 


M'i 

-0. 

-0. 

-0. 

-0. 

-0. 

-0. 

-0. 

-0. 

-0. 

-0. 

-0. 


PENAI   iv 
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/»  I     < '  i  '   I 
A2    I  '      I 

ci  :  i 


0.23797E-01  -0. 
0.2ei61E-01  -0. 
0  . 3  0  S  9  S  E-  0 1         -  0  . 


-0. 
-0. 

-0 


0.10000E  0  1 
0.11  II  OL 
0  .  1 0  0  0  0  E    0  1 


i .  o  o  o  o  o  r  -  o  'j 

L.OOOOOE  -0  5 

1.00  0  0  0003 
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| 


;  c  i  n  r 
u  l 

.' 

Vv3 
PI 

R2 

1 
MS 

I  '. 

rsFii 

III  u 


i  V    MJf'I'i  R 


I  I  PU  i 


LI  (.P   1 


1 


NUK8FR   '5 

UN MARY 

2    LOOP 


12/12/66 

OF  INDIVIDUAL 
5    LOOP   4 


I  OPS 
I  OP 


LOOP 


1  nop 


LOOP   8 


! 


I  I  s 

.02, 
21 
907 

4  2  4  0  0 


02.  J2  0M 
2829..8S0 

79  )  ,  .07.  790 
0  0  0  3  9  8  0  0  .  0  0  0 
000424  00.000 

0  C  0  4  6  0  0  0  .  C  0  0 


ENCEM    V/M-  IABLES 
r  C  T  A I     i  i  S  T 
UC  I  ICN 

I  I  I  O 
H 
l  SI    '' 

I  ,l  II 

!(.  H 
li 
[i    I C 

\  c 

AIM  A      l 

i    2 
tREA    3 

i  I     h 

I      I         i 
i 
il 

R2    CCS1 

i   ■   I 

i  i  •   | 

A 0    CCS  1 

'  I 

■   I 

i  !  1 


5  J  I 

0 

2  74 


000 
064 

400 
000 

000 


511 

0 

274 

25  fl 

84 


ooo;<i 

068 
390M 

010M 


2602.320 

I 
2907.79  0 

3  9)100.000 
42  4  00.000 

4  6  0  00.000 
531.000 M 

0.06  4 
2 74. 390L 
2  5  0.02  0M 
84.990U 


260  2. 
2829. 
29  0  7 


32  0 
860 
790 

00  0 


4  2  4  0  0.000 

4  >  0  0  0.  00  0    460  00 


31  1 

0 

274 

250 

84 


0  0  0M 
06  8 
390L 
0  3  0  M 
99  0U 


2602.320       2602.320 

29.860       2829.86  0 

290  7.79  0       2907.790 

39  30  0.00  0  39800.000 

424  00.00  0  4240  0.00  0 

0  00  46  0  0  0.000 

531.0  00M   51J..OO0M 

0.063      0.068 
274.3901.   274.39  0L 
25  0.040M   25Q.Q5  0M 
84.99  0U    3  4.99  0U 


260  2. 

1  (29, 
2907. 

'00. 
4  2  4  0  0, 
46  0  00. 
521 

0 


3  20 
8  4  0 
790 

00  0 
000 
0  10 
OOOM 
0  6  3M 


274.390L 

2  5  0.0 5 5M 
8  4.99  0U 


2  6  0  2. 

290  7. 
980  Q. 

24  00. 

6000. 

521  . 

0. 

2  74. 

2  5  0  . 

84. 


320 

I 
79  0 

000 
0  00 
000 
000 
Q67M 
390L 
06  OM 
99  OU 


'  129. 

290  7. 
398QJ).. 
42400 . 

4  6  0  0  Q . 

521. 

274. 
250. 
84. 


0.287     0.287 
8340.000  8339. 990U 

I  1405. 5L71 /C66.624 


2  04./  I 
1  5 4.41  B 
103.854 

20  1. 

15  2.014 

1  0  L  .  5  5  2 

0  ,49E     06 

'  .  K.I 

1174.1 
129  J. 249 

176  '  . 
0 
D 
i) 
0 
0 
0 

I) 

D 
0 
0 


0  12 

1  I  J 
0  0  7 

001 

I 

02  3 
0  |  I 


203.  168 
132.  122 
10  0.352 
2^0.813 
14  5.6a 
9 7.^14 

0  .  4  7  E     06 

12.        Q 

1  0  2  6  . 
12  2  6., 
1383.471 

15  149. 8 5  1 
0.030 
0 
0 
0 
0 
0 
0 
0 
0 
0 


128 
00  7 
003 
00  1 
0  0  9 
0  )  I 
02  5 
0  25 
0  3  ) 


0. 

833  9 . 

I  8 5  3  9. 
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18114. 0  9  9 


2602.320 
'   I29.J36.Q 

2907.790 

3980  0.00  0 

4  2  4  0  0 

46000 

521 

0 


260  2. 32  0 
2  829.86  0 
29  0  7.79  0 
39  8  00.00  0 
42400.000 
460  0  0.00  0 
521.000 
0.0  65M 
274.390L   274.390L 
2  5  0.0  9  6M   250.0  9  1A 
84.990U    84.990U 


000 
0  Q  0 
00  0 
0  6  4M 


20  4. 
15  4. 
103. 
201. 
151. 
100, 
0.4  8E 


121 

021 
20  0 
779 
60  8 
8  7 2 
0  6 


8  5  13 
0, 
0, 
0, 
0, 
0, 
0, 
0. 

I), 

0. 
Q, 


0  74 
032 
130 

007 
0  0  3 
001 
0  05 
001 
02.5 
0  2  9 
032 


3  3.  32  2 
10  17.749 
1.2  0  7.69  9 
1  3 3  5." 6 8 
9,5  M,  771 


0.237 
8  339.970U 

13  3  2  3,023 

204.290 

154.  3  76 

10  3.  7  39 

201.951 

15  1.949 

10  1.430 

0.48F    06 

3  3.  34  1 

1022.628 

L  2 12 . 2  4  3 

I  336.611 

85  1  L  3.3  3 


0.287 
8339. 970U 

1811 4.099 

2  0  4.122 

154.02  3 

103,20  I 

2  01.780 

15  1.610 

10  0.373 

0.48F    0  6 

33  ._32_4 

10  1 7.74  9 

I      ;  7,69.9 

I  I  S5.868 

952.1.878 


0 
0 
0. 

0. 

o. 
o, 

0, 
0. 
0, 
0. 


0  3  2 

1  LO 
0  0  7 
00  J 

00  1 
006 
001 
Q  2  4 

0  29 
0  32 


0, 

o 

0, 
0, 
0. 
0. 
0, 


0  32 

1.  3  0 
00  7 

0  0  3 
0  0  1 
0  )5 
0  0  L 


0. 

0 

0. 

0, 

0. 


o  ,.q  L5 
0  .  J  19 

0.O32 


0  '  2 

1  3  1 
0  0  I 
00  s 
0  0  1 

0.00  6 
001 

1 1 .-, 

0.029 

0  .  I       ' 


r'l.l.  I  I     .      S.P.       *       I  0  0  • 

II      II     FACTUM 

(111     Ul     ■•  I  /i  I  I  I"    i  ,v    lull'       25 


-  0 .0  03 

-0.003 

-0.00  I 

-0.003 
0.  062 

-  0  .  0  0  1 
0.  062 

-0.003 

-0.001 

0.125 

0  .  L  2  5 

0.  062 

0.062 

0  .06  2 

82 


PUN    MM  1  R        1 

X 

I 

1 

X 

! 

\\? 

X 

' 

' 

X 

i 

1 

X 

5 

i  ? 

X 

6 

K3 

X 

I 

■ 

X 

! 

CSI    !G 

X 

9 

1  !: 

X 

10 

rSFl  I 

X 

1  L 

TLF.3I 

Y 

1 

ii   l  0L    CQS1 

Y 

2 

PRCCLCTICN 

Y 

J 

FEEC 

Y 

4 

rsFii 

Y 

5 

i  si  2a 

Y 

6 

rSFBC 

Y 

I 

II  11 

Y 

8 

rc2o 

Y 

9 

re  3c 

Y 

10 

c 

Y 

L  I 

P.E/»    0 

Y 

1  2 

a  i-l;    i 

Y 

i 

A     2 

Y 

1'. 

■    ! 

Y 

i  5 

MM.     V, 

Y 

l     (i,l 

Y 

1  / 

l 

Y 

18 

H    c  i  *  i 

Y 

cesi 

Y 

20 

!    CCST 

Y 

M-    i  i    v  | 

Y 

22 

a  o   1 1  s  r 

Y 

23 

M    CO  SI 

Y 

2  4 

A/    CCS1 

Y 

CCS] 

CYDX 

Ml 

l  T3LEN  NUM 


LOOP  MfMLR  25 


l,>/12/66 


INPUT 

2  6  2  9 
2  9  0  / 
391 

4  24  00 
46  0  00 

0 
274 

250 
8_5 

0 

9 

L84  05 

2  0  4 
154 

103 

152 

10  1 
0.49 

33 

1  174 
L29  3 
B763 

0 
0 
0 
0 
0 
0 
0 
0 
0 
0 


VALU!. 
.  MOO 
.8900 
.  /MM 
.  0000 
.0000 
.0000 
.0000 
.0640 
.4000 
.  0  COO 
.0000 

.2f67 
.9999 

.5171 
.  2  6  0  5 
.4179 
.  8  5  3  9 
21 
.0142 
.  5  5  2  2 
3E  06 
.7668 
.  /676 
.8289 
.2466 
.8059 
.0325 
,132  7 
.  0 1)  t  6 

.COM 
.0013 

•J)  0:5  2 

.  C01  3 
18 

.0282 

.  C  J  I  ) 


I  INAL 

26  0  2 

282  9 

2  9  0  / 

398  0  0 

4  2  4  0  0 

4  6  0  0  0 

523 

0 

214 

2  50 

B4 

0 

8  3  3  9 

L  8 1 1  4 

2  04 

154 

103 

201 

151 

100 

0  .  4  8 

33 

1017 

1207 

L335 

5  5  2  1 

0 

0 

0 

Q 

0 

D 

0 

0 

0 

0 


VAI 

.3200 
.  16  0  0 
.  /  9  0  0 
.0000 
.0000 
.0000 
.0000 
.  C649 
.39  0  0 
.0968 
.9S.Q.Q 

.2  ('66 
.9700 
.0094 
.  1219 
.0226 
.2011 
.  n  0  1 
.60  56 
.  8  1 2  8 
4FJ  C6 
,323_9 
.7491 
.6591 
.  577 
.  8  /  8  2 
.0320 
.  1  3  0  2 
.  0066 

I  5 
.  0  0  I  3 
.005/ 
.0013 
.024  4 

89 
.0320 


INPUT 


F  INAL- 
0. 
-0. 0299 
0. 

0. 
0. 
0. 
-10 
0 
-0. 
0 


MGDI I 


S  .  P . 


0000 
0  0  09 
0100 
0968 
0.0100 


-0 

-0 

-291 

-0 

-o 
-o 

-0 
-0 
-0 
-928  5 
-0 
24 

4  2 
13  ) 

-0 

-I) 
0 
0 

-0 
0 

-0 
0 
0 
0 


.0001 
.  02  99 
.41/7 
.13  8  7 
.  39  53 
.65  29 
.14  2  0 
.4046 
.6794 
.4531 
.44  2° 
.9815 
.  8  7  0  1 
.6190 
.072  3 
.0005 
.0025 
.0000 
.0000 
.0000 
.  00  05 
.  0  0  0  0 
.  0  0  06 
.  0  0  0  8 
.0010 


0.0000 

0. 

0.0000 

0. 

0. 

0. 

0, 

0, 
-0. 

0 
-0 


0000 
0350 
0001 
0006 
0009 


-0.0000 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 
0. 


MINI.    LIMIT 
-0. 
-0. 

-0. 

-0. 

Ml. 

-0. 

-0. 

-0. 
>  /4 
24  9, 

84 


3900 
9900 
99.0.Q 


8  3  3  9.9700 
-0. 
-0. 
-0. 

15.0  010 
mO. 
-0. 
-0. 
-0. 
-0. 
-0. 
M). 
-0. 
-0. 
-0. 
-0. 
-0. 
-0. 
-0. 
-0. 
-0. 
-0. 

-o. 

-0. 


' . 

Q  Q 

10  0  0 

90000 
90  0  00 

90'       ' 

1 

279 
255 

85 


1  : 

0  0  0  0 

0  0  0  0 

I 

00  0  0 

00  0  0 
4  100 
0100 
0  10  0 


8  J  39.9 

-0. 
2  79.40  0  0 
-0. 
-0. 
-0, 


■0. 

-0. 
■0. 

-o. 

-0. 

-o. 

-0. 

M 

-o. 
-0. 
-0. 

-o. 

-M 
1  . 


1      LI      [  1 
I   \  i  n  o 

12.50'  00  0  00 
]/.'■■  10  ii  0  00 
5  7.50000000 
w .  ooo 

j  /.    oo  ioooo 

1  .  ! 

0.001 

0, 

0.0  0  062' 

0. 


83 


Y 

1 

Y 

2 

Y 

3 

Y 

4 

Y 

5 

Y 

b 

Y 

7 

Y 

8 

Y 

9 

Y 

1J 

Y 

LI 

Y 

12 

Y 

13 

Y 

14 

Y 

15 

Y 

16 

Y 

1  7 

Y 

18 

Y 

19 

Y 

20 

Y 

21 

Y 

22 

Y 

Pi 

Y 

2'. 

Y 

2 'j 

•0.309 
0  .  1  0  0 
0  .  2 1 S 
•0.162 
•  0  .  1  b  1 
■0.141 


I  6  2 
L52 
14.1 


LS4 
5  39 
707 
0  .627 
0.181 
0  .3o7 
G. 
0.15  5 
0.106 
0  .4  56 
0.108 
0  .  6  9  2 
C  .129 
0  .  1  6  9 
0.150 


21  i1  > 
0E  0  1 
?  E  0  I 
4  E-  0 1 

S,E-T  .0 1 

'  i  -01 
9  £-01 
2E-01 

311-0  1 

3.E--.Q3 

6E    00 

or- oi 

5E-01 
5C  0  2 
fit.- 0  5 

7E-05 
QE-.0.5 
6E-G6 
££=•.0.4 
9E-08 
a.F-04 
5E-05 
4  r  -  o  5 


314 

100 
2  1 9 
17  7 
143 
133 
L79 
144 
1 3  3 


i  a  _ 

514 
0.606 
0.63  0 
0.170 
0  .  3  &  7 
0. 
0.  16  9 

?_4  / 


,4  3  5 


101 
692 
123 
0.14  5 

0.151 


7E-04 

or  oi 
7  hi  oi 
IE- 02 
SC-01 
2E-01 
71.-02 
IE- 01 

3£=_Q1 

3E-0  3 
OE-01 
OE  0  0 
8E-01 
OE  0  2 
8E-0  5 

8E-06 

:i!>06 
IE- 0  6 
7E-04 
9E-0R 
2E-05 
3E-04 
2F-0  5 


0' .  3  1 
0.  10 
0.21 
0.17 
0.3  7 
0.  12 


17 
37 
12 

18 
51 
4  0 
65 
15 
38 


3.3  Err  0  4 
OOE     01 

57F  01 
7  1 E  -  0  2 

0  5F-02 
2  2E-01 
9  7E-02 
75E-02 
2 OE-01 


_03 
01 


4  3F- 
4  0E- 
06E-01 
2PE  0  0 
41£  02 
78E-05 


0.944 

0. 
0. 

0.80  6 
0.8  06 
0.8  06 
805 
806 
306 


0 

0 

0 

0. 

0.83  9 
258 
17  0 
101 
L15 


0. 
■0.  16 


3  4 
32 
9.2 
6  9 
12 
11 
15 


98E-06 
7  6E-0  7 
36E-06 
2  IE- 05 
29E-08 
32C-05 
52E-05 
6riE-04 


0. 
■0 
•0 
•0, 

0. 

0. 


0  8  1 
505 

2  36 
690 


0 
0 

0 
•0 

0.315 

0.62  0 

■0.411 

■0.24  0 


4E-07 


3E-03 
3E-03 
3E-03 
8E-03 
3E-03 
3E-03 

2E--  0-4 

8E-0  1 
4E-04 
7E-04 
4E    01 


3E-07 
4E-07 
0  E  -  0  7 
6E-06 
5E-0  8 
3E-06 
3E-09 
6.E.-09 


■0.4098E-07 
0. 

0. 

0. 

0.83  0  4F-0  3 

0.880  4E-0  3 

0. 


879 

88  0 


302 
2  00 
126 


0 

0 

0 

0 

0 

0 
■0 
■0 

0. 

0. 

0. 

0.438 

0.257 
-0.  754 

0. 

0. 

0, 
-0 


8E-03 
4E-0  3 


8F-01 
9E-04 
OE    0  1 


6E-07 
7E-07 

1E-06 


725 
431 


9E 
2E 


-06 
•0  9 


-0, 

o, 

0 

0, 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 
-0 

0 

0 

0 

0 

0 
-0 

0 

0 

0 

0 


L919E-0  7 


2415E-0  5 


8820 E- 03 


881.1E-0  3 


316  2E-01 

1262E     01 


1652E-07 
7  5  5  2  E-  0  6 


75 78 E- 06 


0 

0 

0 

0 

0 

0 

0 

0 

0 

0 

0 
■0 
•0 
•0 

0.4  95 

0. 

0  .  2  5  0 

0. 

0.. 

0. 

0  .296 

0.205 
-0.785 
-0.971 
-0.108 


0.3581! 

0. 


923 
5  46 
327 
4  05 


v 


52 


7£ 

6E 
2E 
5E 
3E 


0 


15E 


03 
•01 
01 
01 
01 
02 


0E-03 


4E-0  4 
4E-05 
2E-04 
2E-0  4 
5E-03 


■0.2754E 
■0.5 76  IE 
■0.8773E 
■0.279  31." 
■0.5869E 
-.Q._9.080E 
0. 

_o .  ' '65r. 

-0.7995E 

-0.74  76E 
-.0  .  138  IE 

0  .  1 6  1  1.  E 
-  0  .  6  0  2  8  E 

0. 

0.264  OF 
-0  .54  051 
-0.1279E 

0.96 
-0  .  1077F 
_tt*_L9J 
-0  .  179  2F- 
:.Q..3  31  IE 


01  -0 
0 
0  6 
0  3 
0-3 
0  3 


■04 


0  3 
03 

0  3 

07 
04 
0  4 
04 

0? 
0  0 

01 

02 

01 

00 

0  2 

0  Q       0 

0  0       0 

0  l      o 


o. 

0 

0 

o. 

0 

0. 

0, 

-o. 

-0 
0 
0 
0, 

-0 
0 
0 
0, 
0 

o. 

__0 

I 


371  IF 
L196E 
15  72F 
15  26E 

I 


03 
01 
01 
01 
01 
02 


.  L193E-04 
.  449 8 E- 04 

.37 '.'•!' -iV, 

-04 

.4  33  IE- 0  4 


l  Y  C  X     V  J!  T  R  I  X 


Y 
Y 
Y 
Y 
Y 
Y 
Y 
Y 
Y 
Y 
Y 
Y 
Y 


1 
2 
3 
4 
5 
6 
7 

9 
10 
11 
12 
13 


X     1C 

0  .     I C  S  E  -  0  3 

0. 

0. 

0  .  1 0 0 1 E  0  1 

0.10  0  IE  0  1 

OIF  01 

n  o  l  r  oi 

l  o  o  1  E  01 

loo  i r  oi 


X     11 


L  1 7  3  E    0  1 


- C  .  5  18  6  E     0  0 


■0 
0. 
0. 
0. 
0. 
0. 
0. 
0. 

). 

0. 
0. 
0. 
0. 


8568E-03 


84 


I  •'< 
L5 
16 

1  7 
18 
L9 

2  0 
21 
22 
2  1 
24 
25 


0.61 C 4 f     0 0 
-  0  . 1  4  3  2  £    0  4 

0. 

Q.  _ 

0 

0 

0 
-0. 

0 

0 
-0 

0 


1  k  3  2  E     0  4 


106  4E-0 3 

i  ■  '  l    -  0  4 

5 4E-04 

03 

,441  IE- 0  4 


L28  1E-04 
1444E-.Q4 


0 

(1 

0. 

0. 

0. 

0. 

0. 

ior-oi 
o. 
o. 
o. 

0. 


SECCNCS     h 


RUN     I  t  ■  ;   '  hi 


P       illM    NUMBER        1 


.001 


INPUT 


LOOP       1 


[NCEPI  .\Li  \  1     \!£H 

1 

1 

R2 

R3 

MS 

CSI 

TS 

T  Sh  L  I 

l"LF  II 

DE  PEN  CENT    VAF  [A 

TCTAL     Cl.Sf 

0DUCT1GN 

FEED 

rsF  io 

TSF2C 

TSI     1 

TC1C 

rC2G 

re  ii 

Q 

A     0 

A     1 

A    2 

A    2 

CCCL    k 

F    Ci       I 

CCST 

Rl    CCS1 

R2    CCST 

R  3    i  :  -  r 

Ca    CCS] 

AO    CCST 

A 1    ■  i   ■   ' 

A2    CCST 

CCST 

LA8LES 


I     I 


JUM  3ER       9 

SUMMARY 
2    LOOP 


12/ 12/66 

OF     INDIVIDUAL 
3  LOOP      4 


LI  CPS 

I  (    IP 


I  COP 


LOOP        / 


LOOP 


! 


85 


LOOP       9 


!LES 


?.  2  3  1 
1 
368  2 


«/. 


^  3 


36352 
32A39, 

455 
0. 
2  7  4 

2  5  0  , 
5 


U 4    2231.  164M 
2  16    2426.276 
5  59    3682.55  0 
5 2 0 34123 . 5 2 0 
70  036  3  5  2 
25JL2SA3J9 


266 
0.5.4 

400 
000 
000 


4/5. 

0. 

2  7  4. 

25  0. 

8  4 . 


700 

25  0 
266M 
0  56 
"i  g  q  w 

0  L.0M 
9  9  0  M 


2  231 

2  4  2  6 

368  2 
3  4  12  3 
3  6  3  5  2 
3  5)4  39.25  0 
495.266M 


16  4 
2  76 
530 
520 
700 


0 
274 
25  0, 
84, 


05  4 
3  9  0U 
010U 
99  OU 


!1. 164 

2426.276 

3682.53  0 

34123.520 

152.700 

39439.250 

515.266M 

0.056 


27  4 
250, 

84, 


390U 
0  1  0  U 
99  OU 


2  2  31 
2426 
3682 
341  2  J 
3  6  3  5  2.700 
39  439.250 
495. 266M 

0.054 
2/4. 390U 
250. 010U 
84.99  0U 


164 
276 

530 
520 


22  J] 

2426 

3  68  2. 

34123, 

3635  2, 

33.43  9. 

515. 

0. 

274, 

2  50. 


164 

2  7  b 
5  3  0 

52  0 
7  0  0 

25  0 

26  6M 
056 
39  OU 
010U 


,  164 
,276 

,  5  3  0 
520 

DO 


34.990U 


2231 

2  4  26 

3  6  8  2 
3  412  3 
3  6  3  5  2 

33439.  2  50 
505.266M 

0.055 
274.390U 
250.010U 
84.990U 


22  3  1.164       22  5  1.164 
2426.276       2426.276 
3682.5  10        1682.530 
34  123.52  0    34123. 
36  352.  700    36352.  MO 
39439.250     ?,  !  5.9.250 

50  5.266  515.266M       5  0  5.2- 

0.053M  0. 

274.390U      274.3 
25  0.010U       25  0,0 i OU 


2  231.164 

!       '  .   '  f ', 
■ 
'3.520 
36352.700 


84.99  0U  84.990U 


0, 
274 


U 
0  I  0  u 

1 


0.293 
.  !  S9.999 
2  3  7  0  3.  15522358 


0  .  2  9  C 
8 3  3  9.  9 9 0 U 


0.28  9 


0 


210. 

16  7. 

102. 

208. 

165. 

10  0. 

.  4  2  L" 

25. 

1  0  4  5  . 

I 

212  2. 
L  0  i  5  , 
0  . 
0  , 
0. 
0 
0 

o 

0 

0 
0 


6  5  2 
L8J 

1 5  2 
379 

5  5  1 
111 

06 

5  50 
■":  ' 
44  4 
04  7 
2  24 

0  4  2 

1  14 
0  0  5 
)0  ! 

0  02 
0  0  1 
I 
0  2  5 
0  25 


057 
7  06 
IS_£ 
179 

426 


0.290 
8339. 970U    8339. 970U 


209 

16  5 

99 

207 

163. 5 4  3 

9  7  .  1 ',  4 

0.44E    0  6 

3  0. 59  6 

9 3  9.1  7  4 

1082.2  16 

1  5  L  8  .  7  5  8 

9  J 45. 82  5 


23  601.99  1 
2  10.592 

16  7.6  5  C 

102.499 

2  0  8.318 

1  6 5  .413 

10  0.51  5 

0.46E     06 

')  .  /<;  / 

9  0  4.333 

103  7. 54  8 

1806.67  8 


0.289 
8  339.970U 


22250, 

209, 

165. 

98, 

207, 


045 
6  34 
64  7 
909 
353 


8  33' 


29  0 
9  7  0  U 


010 
490 


0 


16  3. 391 

9  6.869 

.48E    06 

32.  /,;l 


0.051 


0 

0 

0 

0, 

0 

o, 

0, 

0. 

0, 

0, 


C  15 
119 

005 
0  0  3 
0  0  2 
0  0  6 
00  1 
02  3 
0  26 
0  4  6 


3613 

0. 

Q 

0. 

o, 

0. 

0, 

0. 

0, 

0. 

0, 


591 

04  2 
124 
0  05 
003 
002 
00  2 
001 

'),'; 

025 
04  3 


821, 

94  4, 

1652. 32  5 

1  3  34  U  7 3 4 

0, 

0 

0, 

0. 

o. 

0. 
0. 
0. 
0. 
0. 


008 
0  1  7 


,  0  39 
.  129 
,00  5 
0  0  3 
002 
008 
0  0  1 
02  0 
02  5 
04  0 


2  3  4  5  6 
210 
167 

102.  117 
2  0  3.216 
165.103 
10  0.127 
0.46E    06 
3  1.78  6 
901.  125 
10  34 
1803 
43  6.9 
0 
0 
0 
0 


2  2  14  4.667 
209.55  0 
1  '-  5  .  4  7  1 
9  8.595 
20  7.  26  8 
16  3.214 
96.5 5 0 
0  .  4  8  E    0  6 
.771 
8  1  .-1  .741 
94  1.678 
1649. 88  1 
14  28  2 
0 


.226 
.  0  34 
.  166 
,  041 
,124 
,005 
003 
0  0  2 
0  0  J 
001 
0  2  2 
0  25 
043 


6  1 6 
0  39 
12.9 

;)  0  5 
10  3 
00  2 
0  09 
0  0  I 
020 
0  .02  1 
<t0 


0.2  39 

3339. 97 OU 

2  3  3.02.A28 

1 1 0 .  3  3  2 

16.7..  212 

10  1.713 

2 0  8.10/ 

164.9/0 

99.  71  / 

0  .  4  7  E    0  6 

3  2. 323 

368. 02 0 

'!'<  /,  110 

1735.23 6 

5  753.  305 

0.041 

!  ».6 


0, 
83 
24844, 
211. 

1  6 ! !  . 
105. 

20O. 

16  7. 

103.. 

0.47E 


2^0 
9  7  0U 
60  3 
439 
421 
672 
171 
193 
7  3  1 
06 


32. 
I  .  / 1  7 
!  '  > 
1  <  2 
854 
04  4 


LO  3  ] 

1773 
- 1  37  2 

1 


o 

i  3 

, ,  4  .  2  !   3 
210. 

I      '.155 
10  1.612 
2.0     .   '    Q 
164.9  I  3 
9  9  ._6 1  5 
0.48E     06 
3.2  • 
841  .107 
, 
I  •'  1  ■■  ,   ■ 

0. 


9  /OU 
4  4  4 


24794. t 
211.405 

.351 

.546 

209. 13« 

16  7.137 

103.  (  :1  'f 

0.47E    0  6 

l 

)  I  il 

!    :  .  0 

■117  I 

1  . 


l 
0 

0, 

I) 

o. 

0. 

0 

0, 

0. 


00  5 
0  )  J 
0  )2 
0  )  5 
0  )  1 
021 
0  24 
0  h  2 


0 
0 
0, 

o 
I 
I 

0. 


'         '  0.    ! 


I    I   I 

0  0  ! 
)  »2 
02  5 


0.0 


I 

o, 
0. 
0, 
0. 

o. 


0  1  ! 


) 
0, 

1 

0. 
0. 

I , 


! 


86 


uuk    NUMBER       2                 I           LEN    NUMBER        1  LOOP    NUMBER     LB                  12/12/66 

RUIV    fvLMMv       *  SUMMARY    OF     INDIVIDUAL    LOOPS 

INPUT          M:no     l0  LOOP     LI          LOOP    12          LOOP    13          LOOP     14          LOOP     15          LOOP     16         LOOP     17          LOl         !' 

INCEPENC1KT.-V./UUAL1LFS                      ■       ,        '        -  ^       2231.164       2231.164       2231.164       2231.164       2231.164       223l7l64       2?3l7lM 

2426    274:242^276  2426.276       .-,26.276      2.426.276       2426.276      2426.27.6      2426,276 

3682J55'3    -682.530  682.530        J682.530       3682.310       3682.3  50       3682.530       3682.530                          0        >682. 
;                                                d]       %2Q    U23    520    34123.520    34123.520    34123.520    34123.520    34123.520    34123.520    34123.520     ^123.520 

*\                                                       5^'70036352    700  6352.700    36332.700    36332.700    36352.700     36332.700    36332.700    36352.700     '6352.700 

p'                                     M43S.25  0    39439.230    3-e39.230     3.9.435*. 25.0     39439.230     39_43Su25il      I  I        '439.2 

1*^6515  lltl*  510.266M       515.266M       510.266M       515.266M       510.266M       512.766M        il2.766          510.2; 

0    0                  0    033  0.053               0.055               0.053               0.055               0.053               0.054  1             0.055M 

274    390U  274.390U       274.390U       274.390U       274.390U       274.390U       274.390U       274.390U 

250.010U  250.010U       250.010U       250.010U       250.010U       25Q.0.10U       250.010U       25Q...01.0U 

5    000                .590U  84.990U          84.990U          84.990U          84.990U          84.990U          84.990U             >.990U         84.990U 


K3 

IA  S  i  55  .  266 

I      ,  il 

,S  274.4  00 

,  |  i  25il._Q0.Q 

i 


,:,:'         s"                                   0    25,3             0.285  0.290  0.289  0.290  0.289  0.290  0.289                CL.2S.9                0.2 

'      '       c  833^999    833^9  7  CU    8339.970U    3339.970U    8339. 970U    8339.970U    8339.970U    8339. 970U    8339. 970U    8339. 97  0U 

l\        1U,LN                          ,7    ,23    22    663    24/^9.554    23174.004  24675.825    23121,176  24605.854    23546.623    2Z6JL6. 92J 

710    652       210    '26  211.369  210.292  211.326  210.254  211.779  210.553          209.894                  .553 

'                                                      '78                  !o  \\\.lll  167.022  169.185  166.-44  160. OH  1*7.570          I                                  '•;/ 

]\\    \                                            10  2    752       1 Clio  2  08.409  101.373  105.249  101.233  105.073  102.354            99.882          102.354 

:i  2oJ:c3o  111:111  W.o^  ^.^  207.973  209.010  20a. 279   ^...^       ^ 

,/,:,  s',1    \tU     832  167.050  164.779  166.960  164.700  166.861  165.331    163.939    163.331 

;                     :      I^Jj  \^\lll  59.373  103.302  99.230  103.1,4  100.368     <      >    10. 

0    47E    06    0.48F     06  0.47E    06  0.48F06  0.47E06  0.48H06  0.47E06  0.48E06       0.48E     06       0.47E06 

fA     n                                               23     5CQ          32.863  32.712  32.853  32,707  32,354  32.702  32.7'>3                         L6                 /,' 

104.;:     ,          840!270  K85.101  839.290  883.760  833.227  882.289  853.464 

ARE*                                           1202.444       363.537  1012.147  362.923  1010.752  061.872  1009.224  97&.4S7                                   * 

AREA    -                                        2    22.0971673.786  1740.128  1672.670  1738,610  1671.457  1736.737  I            •         '          6 

?CCL    J                                      145.224    6767.176  -736.742  7053.236  -482.606  7367.931  199.560  ^32.65110493.423 

CCS                                                         >42             0.041  0.044  0.041  0.044  0.041  0.043  0.0 

0.11',             0.175  0.128  0,179  0,173  0.129  0.12.8  0  .J 

;                                                  0.005             0,003  0.005  0,003  0.005  0.005  0.005 

CCS1                                           0.003            0,003  0.003  0.003  0.003  0.00  5  0.003  0.003 

LST                                             o.oi                 0,007  0.002  ).002  0.007  0,00,  0,002  0.002 

,      |                                                      |             ,.  -O.00Q  o,o,4  -0.000  0.004 

,    CCS!                                             0.                     0.001  '      101  0.001  0.001  0,0)1  0.001 

0.021  0.020  0.1 


A?"  CCS1  0.1]  0.023  0.074  0.07  3  0.024  0.03  3 

A3    CCST  0.051  0.040  0.042  0.040  0.042  0.040 


0.04  7  )41  '.04  0  0.04.1 
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iujn    m; 


IBLEN    NUMBER       1 


!  ! 
I 
.' 

V,  ) 
Rl 
R2 
I 
MS 
CSF 
TS 
TSF 
H.I 

Cf  PE 
I  CI 

!'     [J 

F  E  E 
T  S  F 
TSF 
TSF 
TCI 
TC2 

rc  3 

. 

A  ix  E 
A  R  r 
ARE 
CCC 
F    (. 

i 
R2 

Ch 
AO 
Al 
A  2 
A  3 


INPUT  LC.0 P     1 

P.I  m:i  NT    VJR  [AE1  l  S 

2231.164    2 2 

242  6*276    24 

3'6  82.5  5  9     3  6 

14123.52.0341 

1  6 3  5  2.700363 

3343.9,.  23:iO_3£4 


LOOP 


LOOP 


NUMBER    25 
SUMMARY 

2  0  I. COP 


12/  17/66 
n  INDIVIDUAL 
2  1  LOUP    22 


LOOPS 

I  i  OP 


2  3 


LOOP    24 


I  HOP    25 


)C 

1 1 
i  i 


Nl  ENT    V/K  IABLES 
A  L    C  C  S.I 
C  UC  T  1 0  N 
D 

i  G 

3  C 
C 
D 
C 


A 
A 
A 
A 

L  k 
GST 
■•  i  :  I 
C(  S  I 
C (  ST 
CI  '  I 
'.i  '  I 
COS  3 
CCSJ 

<  i  :  i 

MM 


45  5.266 

0  .  0  5  4 
274.400 

2  b  0  .  0  0  0 

a  '5 .  (i  o  o 


31  . 
26. 
G2. 
23. 
'3  2. 
39. 
10. 
0. 
74. 
50. 
84. 


16  4 
2  76 
53  C 
52  0 

700 
213  0 

2  6  6 

Q  55J4 

3  9  0  U 
0  10U 
9  9  0  U 


2  2  3  1.16  4 

24  2  6.2  76 

3682.530 

34 123. 52  0 

3  6 3 5 2.700 

3.9  4  39...  2  5J. 


510 
0. 
274 
2  5  0 , 
84 


0 

8339 

2  3  7  0  3 

210 

L6  / 

102 

2  0  8 

1 6  5 

1  0  0 

0.42 

29 

L04  5 

1  202 

I      ! 

I  0  4  5 
0 
0 
(] 
0 
0 
0 

I) 

0 
0 
0 


.  2 9  3 

.999 
.11)5 
,i 

.787 
.752 
.  L79 
,551 
.7.71 
E  06 
.  [  '  Q 
.  7  9  3 
.  4  4  4 

/ 
.2  24 
.04  2 
,114 

5 
.003 
.002 

I 
.  I  l 
.5.25 

.  0  5  1 


0. 

8  (39. 

16.. 

209. 

166. 

99. 

2J.7., 

16  3. 

97. 

0.4  7E 

32. 

84  !. 

9  t  (  . 

1689. 

10  3  3  7. 

0. 

Q. 

0. 

0. 

0. 

0. 

0. 

0. 

0. 

0. 


2  39 
970 
9  20 
894 
L9  0 

8  82 
6  15 

9  3  9 

e  5. 9 

06 
54_fl 

1  58 

'  0  3 
18  0 

04  0 
128 
C05 
0  03 
0  0  ? 
0  0  6 
0  0  1 
Q  2  0 
0  2  3 
04  0 


0. 

U    8  339. 

2  3  5  4  6. 

210. 

L6J.. 

102. 

2  0  8. 

16  5. 

10  0. 

0.4  /E 

32. 

86  0. 

I  ?10. 

4  (.  9  8  . 
0. 
0. 
0. 
0. 
0. 
0. 

I). 

0. 
0. 

0. 


266 

0  5  4M 
3  9  0  U 
010U 
9  9  0  U 


2  5 

9  7  0  U 

6  2  8 
553 
570 

3  54 
279 
33  1 

3  6  8 
06 

617 
22  0 
32  0 
69  5 
S  3  7 

04  2 
1  2  8 
0  05 
00  3 
0  0  2 
00  3 
0  0  I 
02.1 
0  2  4 
041 


!231.  164 

2426.2  76 

3  6  8  2. .  5  3  0 

34123.520 

36352.700 

3  14  .3  9.  25  0 

510.266 

0.0  55M 

274. 39 OU 

2  5  0  .  0  1  0  U 

84.990U 


223  1.  164 

2  42  6,2/6 

3  6  8  2.530 
34123.520 
36  35  2 
3943  9 

510 

0 

2  74 


2  231.1  6 '» 
242  6.2/6 
3682.5 30 


70  0 
2  5  Q 
266 
054M 
390U 


250. 010U 
34.990U 


34  12  3, 
3  6  3  5  2 
39439,. 

5  1  0  . 

0, 

2  74, 

250. 

04, 


5  2  0 
7  00 
25  0 

266 
0  5  5M 
39  OIJ 

0. 1 0  u 

990U 


223  1. 164 
2  4  26 . 2  76 

3682.530 
341  23.520 
36352.  !  10 


51  0 

0 

274 

250 

84 


2  2  3  1.16  4 

2426.2  L6 

12.530 

34123.520 

S 6  3  5  2.700 
25.0  I  -  "'  .J23.Q 
266  5  10.266 

Q54M  0.05  5M. 

39 OU  274.39  0U 
0  10U  25Q.0.LO.U 
r1  nu  84  .  99  0U 


0. 

8  3  3  9  . 

2 3  067. 

210. 

166. 

101. 

207. 

16  4. 

99. 

0.4  7E 

32. 

8  5  0  . 

9  7  6. 
1699. 
7  4  0  6  , 

0. 
0. 
0. 
0. 
0. 
0. 
0. 
0, 
0. 
0. 


28.9 

97  OU 
116 

2  16 

3  64 
0  9  0 
93  9 
6  1  9 
08  5 

06 
5  7  7 
35  3 
100 
42  7 
74  6 

04  1 
128 
005 
00  3 
002 
0  0  4 
0  0  1 
02  0 
0  2  3 
041 


0 

8  3  3  9 

2  354  6 

210 

1.6  7 

102 

2  0  8 

16  5 

10  0 

0  . 4  7 

3  2 

i0 

98  6 

L  7 1  0 

4  69 '3 
0 
0 
0 
0 
0 
0 
0 

0 

0 
0 


.2  89 
.97  01. 
.  6  2  9 
.  55  3 
.570 
.354 
.279 
.331 
,  36  8 
t  06 
.  6i  Z 
.22  0 
.  320 
.694 
,912 
.  042 
..12.8 
.005 
.003 
.002 
,  0  0  3 
.00  1 
,  0  :  I 
.  0  2  '* 
.041 


0. 
J    8  3  3  9J . 

23067. 
210. 
16  6. 
10  1. 
207. 
164. 

!9U 

0  ,47E 

52. 

85  )  , 
9  /  6  . 

1699. 

7406. 

0. 
0. 
0. 

0. 
0. 
0. 
0. 

0. 
0. 


289 
970U 

116 

2  16 

8  64 

09  0 
9  39 
6  1 9 
085 

0  6 

5  7  7 

3  3  3 
1  0  0 

4  2  / 
72  8 
041 

1  2.3 
0  05 
00  3 
00  2 
,i  o  4 
001 
0  !0 

02  ! 
041 


0 

3  3  39 

2354  6 

210 

167 

1  0  2 

208 

165 

100 

0.4  7 

3.2 

86  0 

9.8.6 

I  7  1  0 

<t698 

0 

0 

0 

n 
0 

0 
0 

0 

0 


.2  i  ) 
.  970U 
,439 
,55  3 
.5.7.0 
.354 
.2/9 
.33  1 
,  3  6  8 
E  06 
.  6  1  / 
.P.PO 
.320 
.695 
.912 
.0  f2 
.  t  >fl 
.3  15 
.003 
.002 
.00  3 
.  0  0  I 
.021 
.0  24 
.041 


0.  28  ? 
8339. 970U 

2  3  0  6  7,11/ 

2  10.216 

1  i6  ."64 

10  1  .090 

207.939 

164.619 

19  .    I   15 

0  .  4  7  E     0  6 

lZ  / 

85  0.353 

976.  L0.1 


16;)^ 

0, 

0 
) 
0. 

1 

0  . 
0. 


4  2  / 

041 

i  ;  ] 

0  0  5 

I 

o  :o 

'   '  3 
04  I 


I 


I 


88 


I 

I  nop 


NIJMI 

■  I  •   I  \ 
2        I    I 


I 


!        S 


If-     [MO! 

I    ) 


/  ■>  3  /  6  7 

v  rnijAi 

1  »     /i 


■in  >s 


i 


i       •      i 


i 


> 


L'J>J70_._2_72V     1970. 


l  r 


ton    .50  0    •■..■.'  i om 

.      Hi     I ;    '7.  7' 
101.7      l     <i  ''■'.'  70.  78 1 M       016. 

L1iL.-Q0Jlf16^2 .3  .  2  0  0  M  5  ft  3  6  9 .  8  2  0  M5  6  5  5  2  . 
,00.0  0  043610.  0  00  M41,?  8  2.  osp. M3fi96  5. 
'00.  QQQfti     "  '  .    52    £j     S s .  Q •? Q m    -  •. 
•    '   .  ■■  10        '•   .     .  49  3M       47  0. "4 6 5 

n.r  ,>,  ,   >  >.,  Q.120M  0. 

I74.40  0U       2  74.4  0  0U        >  7 '1  . 

?5Q.  OOQ        ?  5  0  .  0  0  0  U  l.O  00  U        "    1  , 


.  '  'i  ') 


)  1    )  u 


-  5  .  0  0  0  U 


8  5. 


'  7 
45 
42 

1 
2  ; 

1  1 
0 

p_q 
00 


'  >490 
2  '  2_0  7 
2M     J  641 

,   ■  1      /  ■ 

7 
7   ! 

OH 

ou 
fiu 


,  195M      !0 
,55      !    40 

.  I   '  1044.9 

.    ''!     ■'■■    >p   ,    I.] 


1 '. 


3  2  .  8 


! 

i  0  0 ! 


8  pj_  \{\  269 
?  0M35  !  7 
)8M6  28  14 


[199 

I 

294M6 


?6R  ' . 

3  3  2  3  . 
2  041. 


0 
2  74 
2.50 


.  0  1  ft 
,_1  0  7  ;-l 
.  4  0  0  U 
.010U 

.  ooou 


i    ' 
'i  1 1  •■• 

'    I  ■ 


162.  0 

n .  0 
/  4 . 4 

5  n  ,  0 

I  5  .  0 


7  J 

8  !  j 
OO'J 
I   )_y 

0  V) 


5  /i.  VI') 

]_,   >79M 

'  74  .  40  01 J 

L5_0_-  LLQU 

15  .OOOU 


0_. 

1  7      . 


; 
1 
I 
1 


238     .        I 

1 2  '1  '   ,  '.         ' 
,      /]    • 
I 
730 . 09 

0.0 
?  74  .400U 
. 
.    I  I  O'J 


'   ■ 


1  >  1  < 

-   ■ 


.  i      !    ' 


>  r :   • 

1      « 

If 

T   I 

t  n 

T  c 

fr  ■  p 

re  ! 

' 

1 
! 

1    \ 

- 

' 

'. 

r 

' 

•    '  '    I  '  ', 


0 .  '.  1 ■', 


D_,    ■:     ■ 
34  0.  !  ■  ■    ■  .    '  I    'i 
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APPENDIX  III 
0ESC3IPTION  OF  CQMPU'IEH  PHQ3IUM  -  riQDEL 

Computer  program  MODEL  may  be  used  for  simulation  studies 
of  the  MEMS  process.   Furthermore,  the  program  may  be  used  to 
obtain  feasible  starting  points  for  SUBROUTINE  MQOEL  which  are 
near  the  boundary  of  the  feasible  constraint  region*.  ■ 
Table  III-l  describes  the  variables  used  in  MODEL  and  the 
SOBrlQUTIJJS  MODEL  used  in  the  optimization  study.  |  a   logic  dia- 
gram for  program  MODEL  is  given  in  Figure  III-l.   Program 
MODEL  is  listed  in  Table  III-2.  * 
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Table  III. 

Explanation  of  ComDuted  Program  Variables  Used 
In  Program  MODEL  and  subroutine;  MODEL 

Symbol 

Explanation 

A 

Constant  used  as  multiplier  to  change  size 

of  independent  variables 

AA 

The  first  term  of  the  Lagnangian  Polynomial 

AS 

The  second  term  of  the  Lagnangian  Polynomial 

AC 

The  third  term  of  the  Lagnangian  Polynomial 

AD 

The  fourth  term  of  the  Lagnangian 

Polynomial 

AL 

Latent  heat  of  vaporization  of  distillate 

water  in  each  effect 

ALS   • 

Latent  heat  of  brine  heater  steam 

ALII 

Flashing  brine  stream  flow  rate  into  effect 

one 

AL2I 

Flashing  brine  stream  flow  rate  into  effect 

two 

AL3I 

Flashing  brine  stream  flow  rate  into 

effect  three 

ALIO 

Flashing  brine  stream  flow  rate  out  of 

effect  one 

AL20 

Flashing  brine  stream  flow  rate  out  of 

effect  two 

AL30 

Flashing  brine  stream  flow  rate  out  of 

effect  three 

AN1 

Number  of  stages  in  effect  one 
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Table  3. 

(Continued ) 

Symbol 

Explanation 

AN2 

Number  of  stages  in  effect  two 

AN3 

Number  of  stages  in  effect  three 

Al 

Average  boiling  point  elevation  in  effect 
one 

A2 

Average  boiling  point  elevation  in  effect 
two 

A3 

Average  boiling  point  elevation  in  effect 
three 

3 

A  constant  used  in  the  equations  which 
calculate  horsepower 

CF 

Salinity  of  feed  stream 

CONST(l) 

through 

co;nst(7) 

Cost  coefficients  used  in  the  terms  of  the 
process  objective  function 

CF 

Heat  capacity  of  flashing  brine  solution 

CSF1I 

Salinity  of  flashing  brine  stream  fed  to 
effect  one 

CSF2I 

Salinity  of  flashing  brine  stream  fed  to 

. 

effect  two 

CSF3I 

Salinity  of  flashing  brine  stream  fed  to 
effect  three 

CSF10 

Salinity  of  flashing  brine  stream  discharged 
from  effect  one 

CSF10 

Salinity  of  flashing  brine  stream  discharged 
from  effect  two 

CSF10 

Salinity  of  flashing  brine  stream  discharged 

from  effect  three 
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Table  3. 

(Continued) 

Symbol 

Explanation 

DTO 

Temperature 

difference  available  for  heat 

transfer  in 

brine 

DTI 

Temperature 

difference  available  for  heat 

transfer  in 

effect  one 

DT2 

Temperature 

difference  available  for  heat 

transfer  in 

effect  two 

DT3 

Temperature 

difference  available  for  heat 

transfer  in 

effect  three 

HP1 

Horsepower  of  the  recycle  pump  of  effect  one 

HP2 

Horsepower  of  the  recycle  pump  of  effect  two 

H?3   • 

Horsepower  of  the  recycle  pump  of  effect 

three 

l.J.K 

Subscripts  and  counters  used  in  program  MODEL 

TLFOI 

Temperature 

indicated  by  (Tj ),  on  Figure  5 

TLF1I 

Temperature 

indicated  by  (Tj )q   on  Figure  5 

TLF21 

Temperature 

indicated  by  (Tj )«  on  Figure  5 

00,  Ul 

Overall  heat  transfer  coefficients  used  in 

the  brine  heater  and  effects  one,  two  and 

U2,  U3 

three  respectively 

XA 

230°  F 

XS 

250°F 

XC 

270°F 

XD 

290°F 

X(l) 

Distillate  produced  in  effect  one 

94 

Table  3. 

(Continued) 

Symbol 

Explanation 

X(2) 

Distillate  produced  in  effect  two 

X(3) 

Distillate  produced  in  effect  three 

X(4) 

Recycle  used  in  effect  one 

X(5) 

Recycle  used  in  effect  two 

X(6) 

Recycle  used  in  effect  three 

X(7) 

Srine  heater  steam  consumption  rate 

X(8) 

Salinity  of  flashing  brine  stream  dis- 

charged from  the  MEMS  plant 

X(9) 

Brine  heater  steam  temperature 

X(10) 

Temperature  of  flashing  brine  stream  fed 

to  effect  one 

X(ll) 

Temperature  of  seawater  fed  to  MEMS  plant 

YA 

Enthalpy  of  saturated  steam  at  temperature 

XA 

YB 

Enthalpy  of  saturated  steam  at  temperature 

XB 

YC 

Enthalpy  of  saturated  steam  at  temperature 

XC 

YD 

Enthalpy  of  saturated  steam  at  temperature 

XD 

Y(l) 

Total  cost  per  1000  gal.  distillate  water 

produced 

Y(2) 

Total  distillate  water  production  rate 

Y(3) 

Seawater  feed  rate 

Y(4) 

_ 

Temperature  of  flashing  brine  stream  leav- 

ing effect  one 
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Table  3. 

(Continued) 

Symbol 

Explanation 

Y(5) 

Temperature  of  flashing  brine  stream 

leaving  effect  two 

Y(6) 

Temperature  of  flashing  brine  stream 

leaving  effect  three 

Y(7) 

Temperature  of  distillate  stream  flowing 

from  effect  one 

Y(8) 

Temperature  of  distillate  stream  flowing 

from  effect  two 

Y(9) 

Temperature  of  distillate  stream  flowing 

from  effect  three 

Y(10) 

Brine  heater  heat  transfer  rate 

Y(U) 

Brine  heater  heat  transfer  area 

Y(12) 

Effect  one  heat  transfer  area 

Y(13) 

Effect  two  heat  transfer  area 

Y(14) 

Effect  three  heat  transfer  area 

Y(15) 

Cooling  water  flow  rate 

Y(16) 

Cost  for  feed  pretreatment  and  pumping 

Y(17) 

Cost  of  brine  heater  steam 

Y(18) 

Cost  of  recycling  brine  in  effect  one 

Y(19) 

Cost  of  recycling  brine  in  effect  two 

Y(20) 

Cost  of  recycling  brine  in  effect  three 

Y(21) 

Cooling  water  cost 

Y(22) 

Cost  of  heat  transfer  area  in  brine  heater 

Y(23) 

Cost  of  heat  transfer  area  in  effect  one 
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T.nble  3.   (Continued) 


Symbol  Explanation 


Y(24)  Cost  of  heat  transfer  area  in  effect  two 

Y(25)  Cost  of  heat  transfer  area  in  effect 

three 
Z  Dummy  variable  used  in  Lagrange  polynamial 
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Table  III-2.   Coraouter  program  MODEL 

HCNSS      JOS 

MCNSS      COMT   8  MINUTES.  9  PAGES. 

MCNSS      ASGN  MJR.12 

HCNSS      ASGN  MGC16 

MCNSS      MODE  GC 

MCNSS      EXEO  FORTRAN. ,.,., , MODEL 

DIMENSICNX(12)  .YI25  1  .CCNSTI7) 

C 

SIMULATION  PROGRAM  MEMS  SEAWATER  DISTILLATION  PLANT 

1  FCRMATi2CX,F 13.6) 

2  FCRMAT(20X,E13.6) 

3  FORMAT (5E13. 6) 

4  FCRMAT16F1G.0) 

READ(1.4) (CONST (I ) .1=1,7) 

6  CONTINUE 

RFADd  .DA 

WRITFI 3,2) A 

K  =  l 

IFiA.NE.l. 1GOTC13 

J=26 

GCTC14 

13  J-l 

14  READ( 1 ,1)  (X(  I ) .1  =  1  .11) 

15  CONTINUE 

*w 

SIMULATION  PROGRAM  MEMS  SEAWATER  DISTILLATION  PLANT 

C 

i  *#*•*** *•»•****#**#*******■*#**#**#*******♦# ******♦***#**•**»** 

c 

CONSTANTS  USED  IN  PERFORMANCE  EQUATIONS 

■  AL=100C. 

8=17. 01723 

ANI=?3. 

AN2=>?3. 

AN3=22« 

CF=0.035 

CP=1. 

XA  =  23C 

X8=250. 

XC=27C. 

XD=290. 

YA=953.8 

YB=945.5 

YC=931.8 

YD=9]7.5 

U0=510. 

Ul=510. 

U2=510. 

U  3  =  5 1  C  . 

C 

ft**************************************** **«#*******#***** 

C 

PERFORMANCE  EOUATIONS 

Y(2)=X(1)+X(2)+X(3) 

Y!3)=Y 
ALlCiY 
Ai_2C  =  A 
A13C=A 
A|_II=Y 
AL2I*Y 
AL3I=Y 

csfic= 

CSF2C= 

CSF1I" 

CSF2I= 

CSF3I 

Y(4)=X 

Y(5)=Y 

Y(6)=Y 

Ai=1.0 

A2=1.0 

A3=0.3 

Y(7)=Y 

Y (8)=Y 

Y (9)=Y 

2=X(9) 

LAGRAN 

AA-YA* 

AS=YB* 

AC=YC* 

AD«YD* 

ALS=AA 

*■#*#** 

YQC}  = 

TLFOI« 

T|_F1I  = 

1/ iCP*{ 
TLF2I= 

1/  [CP*< 
Y(15)= 

1/<C?*( 
DTC=X( 
DT1=X( 
DT2=Yt 
DT3=YC 
Y<11)= 
Y(12)= 
Y<13)= 
Y(14)= 
HP1»X< 

1-EXP1- 
HP2=X! 

1-EXPi- 


f  2)  /(  1 
(3)-X( 
L10-XI 
L2C-X( 
(3) +X( 
(3)+X< 
(3)  +X< 
CF*Y(3 
CF*Y(3 
CSF1C* 
CSF2C* 
X( 8)*( 
( 10 )-( 
!4)-<A 
(5)-< A 
100+ (C 
:;75  +  <C 
201+IC 
( 4 ) -A  1 
(51-A2 
(6) -A3 

GIAN  P 
(Z-XB) 

(2-XA) 
(Z-XA) 
(Z-XA) 
+  A  3t  AC 
***■*-** 

X(7)*A 
XQO)- 
(CP*( < 

V(3)-rX 

<CP*( ( 
YOJ+X 
tYC 10) 

Y(6 )-X 
91-C.5 

ICJ-TL 
41-TLF 
5  1-TLF 
Y( 10)/ 
X! 1 )#A 
X( 2)*A 
X(3)*A 
4)*B*( 
L/(0. 
5)*B*( 
AL/IO. 


-CF/X1S) ) 


X( 
X( 

AL1 

AL2 

L1C 

L2 

30+ 

/CP 

CP) 

CP) 

11+ 

21  + 

31 


X(2) 


+  X<4)  1/AL1I 
+X( 5) 1/AL2I 
X(6) 1/AL3I 
)*ALCG(CSF1C/CSF1I ) 
*ALCG(CSF2C/CSF2I 1 
*ALCG(X(8 1/CSF3I ) 
CSF1C)/(2.*0.0300 ) 
CSF2C>/(2.*0.0347) 
X<8)  1/(2. *0. 0315) 


LYNCMIAL 
*(Z-XC)*(Z-XD) /( (XA-XB)*(XA-XC)*(XA-XD) ) 
*(Z-XC)*(Z-X0) /( (XB-XA)*(X3-XC)*(XB-XD) ) 
*<Z-XS)* (Z-XD) /( (XC-XA)*(XC-XB)*(XC-XD) ) 
*(Z-X3)*  (Z-XC)/(  <XD-XA)*(XD-X3)*(X[>-XC)  ) 
D 

*#**#*#*######*  ######**Hf  •&■«■■«■*#*■*■****■***■*#■***** 
LS 

Y( 10>/(CP*(Y(3)+X(4> ) ) 

AL1C+XI5 ) )*Y(4)+X(1)*Y<7) ) -Y ( 10) ) 

(5)  !  ) 

AL2C+XI6 ) )*Y!5)+tX( D+XI2 ) ) *Y < 8  ) ) -Y ( 10) ) 

!  6  )  )  ) 

+Y!3)*CP*X( 11)-(CP*(AL3C*Y(6)+Y(2)*Y(9) ) ) ) 

(11  )  )  ) 

*(X(10)+TLF0I ) 

F0I-A1-(X(10)-Y(4)  >/<2.*ANl.) 

1I-A2-(Y(4)-Y<  5)1/(2. *AN2) 

2I-A3-(Y(5  1-Y(6) )/(2.*AN3l 

(D70*L>0) 

L/(DT1*U1) 

L/(DT2*U2) 

L/(DT3*U3) 

HXP1-AL/ (0.1104*(X(1C)+460.11) 

1104*(Y( 4J+460. ) ) ) 1 

EXP (-AL/I0. 1104* < Y(4)+460.) ) 1 

1104*(Y( 51+460. 1 1 1 1 
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Sable 

III-2.       (Con't) 

hP3=X(6)*B*(EXP(-AL/(0.1104*(Y<5)+460.))) 

l-EXPl-AL/ (0.1104*1 Yt 6 1+460.  )  )  )  ) 

c 

*#■»#»*»  it*  **********************************  ***************  1 

c 

COST    EQUATIONS 

Y(  16)=C0NST(  1)*Y(3)  /l.E+10 

Yil7)=CONST{2)*X(7)/l.E+10 

Y ( 18 >=CONST( 3 )*HP1/ l.E+10 

Y!  19)=CONST(3)*HP2/1.E+10 

Y(20)=CONST(3)*HP3/1.E+10 

Y(21 )=C0NST(4)*Y( 15) /l.E+10 

Y(22)=CONJST(6)*Y(  11  )  /l.E+10 

Y(23)=CONST(  5)*Y(12) /l.E+10 

Y(2t)=CCNST(5)*Y(13)/l.E+10 

Y(25)=CONST( 5  1*Y( 14) /l.E+10 

Yd >=Y( 16 )+Y( 17)+Y( 18)+Y(19)+Y(20)+Y(21 ) 

1  +  Y ( 22 ) +Y ( 2  3 ) +Y ( 24 ) +Y ( 25 1 +CONST ( 71 /l.E+10 

IF(K-2)91,110,91 

91 

DC93I=1,3 

IF(X( I ) ) 104,104,94 

93 

CONTINUE 

04 

DC96I=1,15     . 

IF! Y( I ) ) 104,104,96 

96 

CONTINUE 

97 

DC99I=1 ,8 

XI  I  )=A*X(  I ) 

99 

CONTINUE 

X  (31=8340. -X(2)-X( 1) 

IFI25-J) 104,102,102 

102 

J=J  +  1 

GOT015 

104 

D0106I=1 ,8 

X ( I )=X( I ) /A 

106 

CONTINUE 

X(3)=8340.-X(2)-X( 1) 

K  =  2 

G0T015 

no 

WRITE (3, 2) (CONST!  I )  ,  1  =  1  ,7) 

WRITE (3,1)  <X( I  1 ,1  =  1 ,11 ) 

WRITEI3.2) (Y( I ) ,1=1,25) 

WRITF(3,3)AL10,AL20,AL30 

WRIT-:(3,3)AL1I  ,AL2I  ,AL3I 

WRITE(3,3)CSF10,CSF20 

WRITE(3,3)CSF1I .CSF2I .CSF3I 

WR:TE(3,3)A1,A2.A3,ALS 

WRITE(3,3)TLF0I ,TLF1 I ,TLF2I 

WPITF(3.3)DTO,OT1,OT2,DT3 

WRITE(3»3)HPl,HP2,HP3 

G0T06 

END 
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Table   III-2.      (Con't) 


MONSS 

EXEQ 
CALL 

LiNKLCAD 
MODEL 

MONSS 

EXEQ 

MODEL. MJB 

DATA 

• 

17650.0 
208000000. 

2500000. 

57590000.0 

1. 
2602.32 
2829.86 

2907.79 

36961. 

39376. 

42719. 

510. 

.065 
274. A 
250. 
85. 

5984.0 

239700.0 

376000.0 

PART    II. 
OPTIMIZATION   OP  A  MULTI-STAG  3   AERATED  LAGOON   BY  THE 
DISCRETE  MAXIMUM  PBINCIPLE 
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1.0  INTRODUCTION 

Pollution  abatement  has  become  a  subject  of  increasing  concern  both 
to  the  technical  expert  in  this  area  and  to  the  common  citizen.   Sources 
of  pollution  range  from  plant  life  to  large  manufacturing  complexes.   The 
tools  for  combatting  pollution  range  from  plant  life  to  man  made  pollution 
control  devices. 

In  most  industries  it  is  not  economically  feasible  to  prevent  waste 
formation.   As  a  result  these  industries  must  concentrate  on  destroying 
or  disposing  of  wastes  once  they  have  been  formed.   The  chemical  process 
industry  and  the  petroleum  industry  serve  as  examples  of  this  type  of 
industry. 

One  means  of  measuring  the  strength  of  a  pollutant  is  to  determine 
the  total  amount  of  oxygen  that  is  required  to  reduce  the  pollutant  to  a 
harmless  state.   Some  pollutants  react  directly  with  oxygen.   The  rates  of 
reaction  of  these  pollutants  are  usually  quite  rapid.   Other  pollutants  are 
degraded  by  bacteria  or  microorganisms  in  an  oxygen  enriched  environment. 
Microorganism  feeding  processes  usually  occur  at  a  considerably  slower 
rate  than  direct  oxidation  reactions. 

This  work  deals  with  the  modeling  and  optimization  of  a  system  for 
treating  petroleum  refinery  waste  water.   The  system  or  process  considered 
is  an  aerated  lagoon.  In  general  the  process  consists  of  introducing 
waste  water  solutions  into  a  large  body  of  water  wherein  they  are  degraded 
sufficiently  to  allow  the  effluent  stream  to  be  discharged  from  the  refinery. 

The  lagoon  model  used  in  this  paper  is  partially  patterned  after  the 
aeration  basins  of  American  Oil  Company's  Sugar  Creek  Refinery  which  is 
located  in  Sugar  Creek  Missouri.   Several  of  the  constants  used  in  the  model 
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were  calculated  from  operating  data  from  this  particular  aerated  lagoon. 
A  detailed  description  of  the  Sugar  Creek  aerated  lagoon,  its  operating 
characteristics,  and  the  nature  of  the  wastes  it  treats  is  given  by 
Stroud,  Sorg,  and  Lamkin  (1)  and  Burkhead  (2).   A  short  description  of 
the  Sugar  Creek  Lagoon  is  given  here.  Waste  water  is  first  introduced  into 
a  pond  which  has  an  oil  skimmer  trough  at  its  outlet  for  removing  any 
surface  oil  slick  which  forms.  The  effluent  from  the  oil  skimming  pond 
is  introduced  by  gravity  flow  evenly  across  the  inlet  of  the  first  aeration 
basin.   The  first  aeration  basin  is  approximately  713  feet  in  length  by  120 
feet  in  width  with  a  depth  of  10  feet. 

The  effluent  from  the  first  basin  is  introduced  by  gravity  flow  into 
the  second  aeration  basin.   This  basin  is  approximately  700  feet  in  length 
by  120  feet  in  width  with  a  depth  of  10  feet.   The  effluent  from  this 
basin  is  introduced  by  gravity  flow  into  the  second  aeration  basin.  This 
basin  is  approximately  700  feet  in  length  by  120  feet  in  width  with  a  depth 
of  10  feet.   The  effluent  from  this  basin  is  introduced  into  a  settling 
basin  before  it  is  discharged  from  the  refinery. 

In  the  first  aeration  basin,  three  mechanical  surface  aerators  which 
are  driven  by  60  horsepower  electric  motors  are  mounted  on  steel  platforms. 
The  platforms  are  positioned  down  the  basin  center  line  which  is  parallel 
to  the  over  all  direction  of  waste  water  flow.   The  second  aeration  basin 
has  three  15  horsepower  surface  aerators  positioned  down  its  center  line. 

Recently  four  more  20  horsepower  aerators  have  been  added  to  the  first 
basin.  These  aerators  are  currently  located  along  a  line  which  is  adjacent 
to  and  parallel  to  the  inlet  baffle.   However,  these  aerators  are  not  on 
immobile  platforms  as  are  the  other  aerators.   The  aerators  oxygenate  and 
keep  the  waste  water  in  the  basins  mixed. 
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The  optimization  goal  is  to  determine  the  size  of  lagoon  and  the  sizes 
and  positions  of  the  aerators  required  to  achieve  a  specified  waste  water 
conversion  such  that  the  total  cost  of  the  aerated  lagoon  is  minimized. 

The  lagoon  is  considered  to  be  a  stagewise  process.   The  optimum 
aerator  horsepower  and  the  lagoon  volume  required  at  each  stage  are  de- 
termined by  a  discrete  version  of  Pontryagin's  maximum  principle  as 
elucidated  by  Fan  and  Wang  (3). 

2.0  DESCRIPTION  OF  THE  LAGOON  MODEL 

In  this  section  the  kinetic,  flow,  and  economic  models  of  the  process 
are  developed. 
2.1   Ideal  Component  Assumption  and  Kinetic  Model. 

In  the  lagoon  model  a  real  waste  solution  which  may  have  many  types 
of  impurities  is  assumed  to  be  composed  of  wastes  which  fit  into  one  of 
three  categories.  The  impurities  in  each  category  are  further  assumed  to 
act  as  a  single  idealized  component. 

The  first  idealized  impurity  is  a  mixture  of  organic  and  inorganic 
impurities  which  can  be  degraded  to  harmless  products  if  it  remains  in  the 
presence  of  degrading  aerobic  bacteria  in  an  oxygenated  environment  for 
a  sufficient  length  of  time.  This  component  is  commonly  measured  in  terms 
of  biological  oxygen  demand  (BOD)  which  is  defined  by  Eckenfelder  and 
O'Connor  (4)  as  "...  that  quantity  of  oxygen  required  during  the  stabilization 
of  decomposable  organic  matter  and  oxidizable  inorganic  matter  by  aerobic 
biological  action." 

The  kinetic  expression  which  will  be  used  to  describe  the  rate  of 
destruction  of  the  BOD  component  will  now  be  derived. 
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Inspection  of  equations  given  by  Grieves,  Milbury,  and  Pipes  (5)  shows 

that  the  rate  of  formation  of  aerobic  microorganisms,  r„ ,  can  be  written  as 

B 

dCB    kMCBx 

'■B   dt     (K  +  x)                                   (1) 

where, 

C  =  aerobic  microorganism  concentration 

M 
k  «  growth  rate  constant 

x  »  BOD  concentration 

K  =  Miachaelis-Menton  constant 

If  it  is  assumed  that  the  organism  population  increase  is  proportional 

to  the  increase  in  BOD  concentration,  that  is 

dC  -  -Y  dx                                          (2) 

or 

^=-Y^                                          (3) 
dt       dt                                          w; 

or 

rB  "  "Y  rx                                          (« 

where 

Y   =  a  yield  constant  (lbs.  microorganisms  formed  per  lb. 

BOD  consumed) 

-r  ■  BOD  reduction  rate 
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the  BOD  reduction  rate  can  be  written  by  combining  Equations  (1)  and 
(4)  as 


M 

. k  s  * 

rx  =  Y(K+x) 


(5) 


To  make  the  eventual  optimization  problem  easier  to  solve,  it  is 
assumed  that  Equation  (5)  may  be  approximated  by  a  pseudo  first-order 

kinetic  expression  in  terms  of  BOD.   The  approximation  requires  that 

M 
the  quotient,  k  C  /(Y(K-toc)),  remain  nearly  constant  throughout  the  portion 

of  the  lagoon  which  decomposes  BOD,  that  is 

-rx  =  k  x  (6) 

where 

M 
1    k  C 
k   " B   ~*  constant 

Y(K-hx) 

The  validity  of  the  assumption  made  to  obtain  Equation  (6)  will  be 
verified  after  the  optimization  study  has  been  completed. 

The  second  idealized  waste  component  is  assumed  to  be  a  mixture  of 
organic  and  inorganic  compounds  which  pass  through  the  lagoon  without 
being  directly  oxidized  or  acted  upon  by  microorganisms.   This  component 
does  have  an  oxygen  demand;  however,  neither  direct  oxidation  nor  attack 
by  microorganisms  while  the  component  is  in  the  lagoon  will  reduce  the  oxygen 
demand.   Consequently,  since  the  component's  oxygen  demand  is  not  reduced 
in  passing  through  the  lagoon,  it  is  called  the  nondegradeable  component. 
The  kinetics  of  decomposition  of  certain  detergents  are  such  that  they  are 
effectively  nondegradeable  in  an  aerated  lagoon. 

The  third  idealized  component  is  assumed  to  be  a  fast  reacting  in- 
organic component  with  the  following  characteristics:   It  is  oxidized 


107 

directly  by  oxygen  in  a  relatively  fast  reaction.   The  reaction  rate  is 
assumed  to  be  so  fast  that  it  is  limited  by  the  rate  of  mass  transfer 
and  the  rate  of  addition  of  oxygen  by  the  oxygenation  equipment.   The 
rate  of  oxidation  is  primarily  limited  by  the  rate  of  diffusion  of 
oxygen  into  the  waste  solution.   The  primary  example  of  this  type  of 
component  in  refinery  waste  water  is  hydrogen  sulfide. 

The  theoretical  oxygen  demand  (TOD)  exhibited  by  the  three  components, 
that  is,  the  sum  of  the  oxygen  demands  of  all  three  components,  is  of 
primary  interest  in  this  work.   The  chemical  oxygen  demand  (COD)  includes 
the  biological  oxygen  demand  and  the  oxygen  demand  of  the  fast  reacting 
inorganic  component. 

Further  description  of  the  lagoon  model  can  be  facilitated  by  referring 
to  Figure  1.   The  values  of  BOD,  COD,  and  TOD  are  shown  at  various  points 
along  the  lagoon  aeration  basin  in  the  figure.   At  the  inlet  of  the  aeration 
section  of  the  lagoon  the  TOD  will  include  the  oxygen  demands  of  the  BOD 
component,  the  nondegradeable  component,  and  the  fast  reacting  inorganic 
component.   At  the  outlet  of  the  inorganic  reduction  section,  the  TOD  will 
be  equal  to  the  sum  of  the  oxygen  demands  of  the  BOD  and  nondegradeable 
components,  that  is,  in  the  lagoon  model  all  of  the  directly  oxidizable 
inorganics  are  removed  in  a  section  at  the  beginning  of  the  lagoon  whose 
volume  should  theoretically  be  determined  by  a  diffusion  limiting  rate 
equation.   In  practice  the  aerators  for  this  section  are  set  as  close  as 
practically  possible  to  the  inlet. 

In  the  lagoon  model  the  inorganic  reduction  section  is  located  at  the 
inlet  because  it  is  desired  to  have  as  large  a  portion  of  the  lagoon  as 
possible  for  BOD  reduction.   This  method  of  aeration  is  chosen  because  it 
is  assumed  that  no  aerobic  bacteria  can  survive  in  the  oxygen  deficient 
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waste  water  which  contains  even  a  small  amount  of  fast  reacting  inorganic 
component. 

At  the  outlet  of  the  aeration  sections  of  the  lagoon,  the  TOD  will 
be  the  sum  of  the  oxygen  demands  of  the  nondegradeable  component  plus  the 
BOD  of  the  outlet  water. 

2.2  BOD  Material  Balance 

The  sequential  arrangement  of  aerators  in  the  Sugar  Creek  lagoon 
suggested  that  the  BOD  reduction  section  of  the  lagoon  could  be  modeled 
by  a  series  of  ideal  backmix  reactors, 

A  steady  state  BOD  material  balance  about  the  nth  ideal  backmix 
(completely  mixed)  reactor  for  the  biodegradable  component  can  be  written 


Qx11"1  -  Qxu  +  r.  Vn  =  0  (7) 


where  Q  is  the  volumetric  feed  rate  and  V   is  the  reactor  or  stage  volume. 
The  volume  of  the  nth  stage  of  the  BOD  reduction  section  can  be  calculated 
by  combining  Equations  (6)  and  (7)  to  obtain 


.,  n-1    a, 

Q(x    -  x  )  ,  . 

,1  n  C8) 

k  x 


2.3  Aerator  Motor  Size  Equation 

It  is  assumed  that  Equation  (8)  is  valid  as  long  as  a  certain  minimum 
oxygen  concentration  is  maintained  in  each  stage.   The  size  of  the  electric 
motor  required  for  maintaining  this  minimum  oxygen  concentration  is  directly 
proportional  to  the  product  of  the  waste  water  flow  rate  and  the  difference 
in  BOD  between  the  stage's  inlet  and  outlet  streams,  that  is, 

.,  n-1    n. 
pn  .  9<*    ■») 


no 


where 

n 
p   «  aerator  electric  motor  size  at  stage  n 

Q   ■  waste  water  flow  rate  through  the  lagoon 

n-1 
x      -  oxygen  demand  of  waste  solution  flowing  into  stage  n 

x  -   oxygen  demand  of  waste  solution  flowing  out  of  stage  n 

Ro  "  oxy§en  transfer  rate  constant 

E   =  mechanical  efficiency  of  aeration  unit 

2.4  Economic  Model 

In  connection  with  designing  an  aerated  lagoon,  the  sizes  and  locations 
of  the  aerators,  the  number  of  aeration  stages  to  be  used,  and  the  lagoon 
volume  should  be  determined  by  economic  considerations. 

The  objective  cost  function  used  in  this  study  takes  into  account  both 
the  initial  equipment  costs  and  the  operating  costs  for  the  life  of  the 
equipment.   A  present  worth  objective  function  of  the  form  used  by  Hwa  (6) 
is  used  to  do  this.  The  initial  costs  considered  are  the  aerator  motor 
costs  as  a  function  of  horsepower,  the  lagoon  volume,  and  the  costs  of  the 
aerators  and  supporting  platforms.   Only  one  operating  cost,  the  electric 
power  cost  for  aerator  operation,  is  considered  to  have  a  significant  effect 
on  selecting  the  equipment  size  and  lagoon  volume. 

Bauman  (7)  has  proposed  the  following  relation  which  can  be  used  to 
obtain  the  initial  cost  of  an  explosion  proof  induction  motor,  (C  n).,  as 
a  function  of  the  motor  size,  Pn, 


(c,a),    =     /S(Pn) 


"I  '1  "  /•"•*  '  (10) 

where 

°^  -  aL   for  l^p'i  20 

^  =  a[   for  20    <  Pn±200 


» 

fi 

«  C,  for 

1  i  Pn  &  20 

ft 

i 
»  C,  for 

20  <■   Pn  £-200 

The  symb 

3 Is  a. 

i 

t 
and  C1  are  constants. 

The 

initi 

ll  cost  of 

the  aeration  basin  is  taken 

as 

the 

sum  of 

land 

real 

estate  cost,  the  cost  of  digging  the 

basin,  and 

the 

cost  of  1 

aying 

rock 

sid 

mg  to 

prevent  erosion.  The  unit 

land,  digging, 

anc 

rock 

lining 

costs 

are  represented  by 

C_,  C  ,  and  C  respectively 

The  aeration  basin 

cost 

for 

each  stage  or  aerator,  (CT);>>  is 

calculated 

by 

Wj*J 

2  -  C2j2w 

+  C3J2wh  +  C4(£w  +  2h^  ) 

(11) 

where 

■ 

St  • 

length  of  aerator's  basin 

w  = 

width  of 

aerator's  basin 

h  = 

depth  of 

aerator's  basin 

To  re 

late 

lagoon  geometry 

to  aerator  requirements,  the  ratio 

of 

J2:w 

:h  =  23: 

12:1 

(12) 

is  us 

ed. 

This 

is  approximately  the  ratio 

of  the  dimensions 

of  the 

Sugar 

Creek 

aeration 

basis. 

Equation 

(12)  y 

ields 

w  - 

11  0 
23 

(12a) 

and 

h  = 

Si 

23 

(12b) 

By  de 

finition 

Vn  = 

iwh 

(12c) 
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Substitution  of 

Equations  (12a)  and  (12b)  into  (12c)  gives 

vn  = 

-,  i3                               (124) 

23" 

Solving  Equation  (12d)  for  J£   gives 

SL3    = 

232  Vn 

12                                  (12e) 

2 

S0lving  Equation  (12e)  for  £  gives 

J22  - 

f232" 

2/3      n  5/<! 

(Vn)2/3                      (12f) 

12 

Substitution  of 

Equations  (12a)  and  (12b)  into  Equation  (11)  gives 

«& 

"   C2  I*'  +  C3  ~2  ^     +     C4  M^    (U) 

23 

Substitution  of 

Equations  (12e)  and  (12f)  into  Equation,' (13)  and 

simplifing  gives 

the  aeration  basin  cost  as  a  function  of  the  basin 

volume 

(cn) 

-     C3Vn  +  C5   (vV2                    (13a) 

where 

a2  " 

2/3 

and 

s  - 

2I  (^]2/3(6C2+  V 

The  aerator 

turbine  and  supporting  stand  is  assumed  to  have  a  constant 

initial  cost,  C, 

irregardless  of  the  aerator  motor  size. 

Since  the  e 

.ectrical  power  cost  is  paid  over  a  period  of  years,  it  is 

necessary  to  estimate  the  present  worth  of  this  money.   The  operating  cost 

(CQ),  is  given  by 
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(C*)   - 

i=M 
C7  T  YZ      (1  +  r)"1 

PU                     (14) 

where 

C-  =  electrical  power  cost  (assumed  constant) 

T  ■  operating  hours  per  year 

M  ■  life  of  system 

r  =  annual  interest  rate 

Equation  (14)  reduces  to 

(C»)  -  C8  Pn                                         (15) 

where 

M 
i=l 

The  total  cost  for  each  stage  of  the  organic  degradation  section  of 

the  aeration  basin  is  then  written  as  the  sum  of  the  costs  given  previous 

ly: 

Gn   -    (c°)1   +   (c»)2   +    (eg)    +  c5 

or 

Gn  =  Cfi  +  C5(Vn)  2  +  C3Vn  +  Cg  Pn  +  /9(Fn)U                      (16) 

The  objective  function  for  minimizing  the  cost  of  the  organic  de- 

gradation section  of  the  aeration  basin  is  then  defined  as 

n=N 

S  "  T.      a"                                        (17) 
n=l 

where 

S  ■  total  cost  of  BOD  degradation  section 

114 


total  number  of  aerators  or  stages 


3.0   PROCESS  OPTIMIZATION 


3.1  Development  of  the  Performance  Equations 

To  restate  our  objective,  the  purpose  of  this  study  is  to  determine 
the  volume  per  stage,  the  aerator  motor  size,  and  the  number  of  stages 
required  for  the  BOD  degradation  section  of  an  aeration  basin  which  is 
used  to  achieve  a  specified  waste  reduction  while  at  the  same  time  mini- 
mizing the  cost  of  the  total  system. 

The  maximum  principle  is  used  to  perform  the  optimization  because 
it  provides  a  systematic  method  for  optimizing  multi-stage  processes. 

In  anticipation  of  the  forms  of  the  equations  that  may  be  used  in 
the  maximum  principle  solution,  some  of  the  previous  equations  will  be 
combined  and  rearranged.   Equation  (8)  can  be  rearranged  to  give 

n-1 
Xl  ,   , 

*1  ■  ? (18) 

1  + vn 

Q 

Substitution  of  Equation  (18)  into  Equation  (9)  and  rearranging  gives 


n-1 
xl 


R  E  (    *    +  J_   )  (19) 

°    k'  V"      Q 


The  electric  motor  size   can  be  eliminated   from  the  stage   cost  equation 
by  substituting  Equation   (19)    into   Equation   (16) 

n-1 

Gn     =     C,   +  C.    (Vn)    2  +  C,   Vn  +  C0  1 

to  j  Jo  — - — ~~    


RuE     t +  ) 

k'V11  Q 
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ft 


oL 


n-l 

*1 


R  E   (- 


k-vn   +   Q    ) 


(20) 


A  new  variable,   x2,   which  is    equal   to   the  sum  of   the   costs   of  all 
the  stages   up    to   and   including  stage  n  is   defined  next  as 


where 


x"-1  +  Gn,  x°  =  0 


(21) 


n-l 


the  sum  of  the  costs  of  all  stages  up  to  stage  n 


G    -  the  cost  of  stage  n 

By  defining  the  stage  volume  of  the  nth  stage  as,  0n,  a  change  in 
notation  to  conform  with  the  notation  given  by  Fan  and  Wang  (3)  can  be 
performed.  It  is  hoped  that  the  change  in  notation  will  make  it  easier 
to  follow  through  the  algorithm  which  is  used  to  solve  this  problem.   In 
terms  of  the  new  notation,  the  dependent  variables  which  are  called  state 
variables  are  denoted  by  the  letter  x.   The  superscript  n  on  a  state 
variable  indicates  that  it  is  the  result  of  a  decision  made  in  the  nth 
stage.   Independent  or  decision  variables  are  denoted  by  the  Greek  letter 
8.   The  superscript  n  on  a  decision  variable  indicates  that  it  is  a 
decision  which  is  made  at  the  nth  stage.   The  BOD  reduction  section  of 
the  aerated  lagoon  process  can  be  visualized  in  terms  of  the  discrete 
maximum  principle  by  referring  to  Figure  2. 

Given  the  values  of  the  state  variables  entering  a  stage  and  the 
values  of  the  decision  variables  at  that  stage,  the  value  of  a  state 
variable  leaving  the  stage  is  calculated  by  using  its  transformation 
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equation  (Equation  (18)  for  x.  and  Equation  (21)  for  x,).   In  general, 
the  transformation  equations  are  the  performance  or  constraint  equations 
of  the  process. 

3.2  Statement  of  the  Optimization  Problem 

The  transformation  and  economic  equations  which  are  needed  for  the 
optimization  study  can  be  rewritten  in  terms  of  the  new  notation  and 
summarized  as  follows.   The  first  number  to  the  right  of  each  equation 
gives  the  original  equation  number  from  which  the  equation  is  obtained. 
The  BOD  material  balance  equation  which  is  the  transformation  equation  is 


n-1 
Xl 


i  +  k  en 


T(x^S  6n) 


(18)  -  (22) 


where  T(x.   ;  8  )  is  the  symbolic  notation  used  for  the  transformation 
equation.   The  cost  at  each  stage  is 


G(x^-1;  9n) 


=  C6  +  C5  (8n)  2  +  C3 


+  C0 


+  $ 


^ 


n-1 
Xl 


R  E  (- 


•-JJ 


1-1 


R  E  ( 


k'8 


-Li 

Q 


(20)  -  (23) 


where 


esi- 

al 

for 

«U 

»i 

for 

ft  = 

ci 

for 

fi- 

< 

for 

1  <  P  <   20  horsepower 

20  <  Pn  <  200  horsepower 

I   <_  P      <  20  horsepower 

20  <  20  _<  200  horsepower 


lis 


and  the  accumulated  cost  is 


x^  =  x*"1  +  G(xf  S  9n),  4     -   0  (21)  -  (24) 


The  objective  function  which  is  to  be  minimized  is 


(17)  -  (25) 


The  transformation  equation,  Equation  (22),  may  be  used  to  calculate 
the  BOD  concentration  at  the  outlet  of  the  nth  stage  given  the  inlet  BOD 
concentration  and  the  stage  volume  (the  decision  variable)  at  that  stage. 
Equation  (23)  gives  the  cost  of  the  nth  stage  in  terms  of  the  inlet  BOD 
concentration  and  the  stage  volume.   Equation  (24)  gives  the  total  cost 
of  all  stages  up  to  and  including  the  nth  stage.  The  total  cost  of  the 
BOD  reduction  section  of  the  lagoon  is  given  by  Equation  (25)  as  the  sum 
of  the  costs  for  each  stage  of  the  BOD  reduction  section. 

Equations  (22)  through  (25)  are  the  performance  equations  for  a  one- 
dimensional  multistage  decision  process  which  is  defined  by  Fan  and  Wang 
(3)  as  a  "  ...  process  which  can  be  completely  characterized  for  the  purpose 
of  optimization  by  a  single-state  variable..."   The  stage  volume  is  the 
decision  or  state  variable  for  the  lagoon  process. 

3.3   Computational  Procedure  -  The  Discrete  Maximum  Principle 

Fan  and  Wang  (3)  have  derived  a  necessary  but  not  sufficient  recur- 
rence relation  which  can  be  used  to  calculate  the  optimal  values  of  state 
and  decision  variables  for  many  one-dimensional  processes. 


-„/  n-1  „n,     n„.   n  .n+1. 
3G(x    ;  e  )     3G(x  ;  6    ) 

„n+lN     „„,  n   n+lN 


3 


,n+l      ST(xj;  e"Ti)     300^;  8"r*) 


„,  n-1  ,  n.     ,„,  n  ,n+l,       .  n  .  n 

r(x   ;  -j   )    3i(x1;  6   )       3x1  3x1 

~  .  „n+l 

36  38 


(26) 
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Equation  (26)  has  the  general  form 

gCx?"1;  9n)   -  f(x";  9n+1)  (26a) 

where 

,  n-1  niiv    ,  -/  n  An+1* 
g(x   ;  9  )  and  ffcCjl  9   ) 

are  the  left  and  right  hand  sides,  respectively,  of  equation  (26). 

One  solution  procedure  for  a  problem  with  fixed  end  points,  i.e., 
known  inlet  and  outlet  BOD  concentrations,  is  to  start  at  the  first 
stage  of  the  process  and  assume  a  value  of  the  decision  variable  for 
this  stage,  8  .  With  Equation  (22),  the  transformation  equation,  x,  is 

calculated.   Next,  Equation  (26a)  is  used  to  calculate  the  optimum  value 

2 
of  the  decision  variable  for  the  next  stage,  8  ,  for  the  assumed  value  of 

8  .   Repeated  application  of  the  transformation  equation  and  the  optimum 

recurrence  equation  for  all  of  the  stages  of  the  process  will  yield  a 

value  of  the  outlet  BOD  concentration,  x,.   If  this  value  of  x^  is  equal 

to  the  required  value,  the  problem  has  been  solved.  If  not,  another 

value  of  8  is  assumed  and  the  above  procedure  is  repeated. 

The  optimization  calculation  was  carried  out  using  the  algorithm 

described  above.   However,  in  carrying  out  the  computations,  the  electric 

motor  size  at  each  stage  of  the  process,  P  ,  was  used  as  the  decision 

variable,  8n,  in  place  of  the  stage  volume  in  order  to  simplify  the 

computations.   The  performance  equations  which  were  used  to  perform  the 

optimization  are  given  in  Appendix  I.   Theoretically,  either  the  set  of 

performance  equations  with  stage  volume  as  the  decision  variable  or  the 

set  of  performance  equations  given  in  Appendix  I  which  have  the  stage 
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electric  motor  size  as  the  decision  variable  will  work  equally  well  when 
performing  the  optimization  calculations.   Computer  program  OPT  which  was 
programmed  to  perform  the  optimization  is  explained  in  Appendix  I. 

4.0  RESULTS  AND  CONCLUSIONS 

The  constants  which  were  employed  in  the  numerical  solution  of  the 
problem  are  given  below 

a.   =  0.53    for   1  <  Pn  <  20    horsepower 
«j  =  1.08    for  20  <     Pn  <  200   horsepower 

a2  =  0.6667 

(BOD).  =  0.0109025  lb.  oxygen/ cu.  ft.  solution 

(BOD)  =  0.0026166  lb.  oxygen/cu.  ft.  solution 

(COD).  "  0.0290941  lb.  oxygen/cu.  ft.  solution 

(COD)  =  0.0090958  lb.  oxygen/cu.  ft.  solution 

Subscripts : 

i  =  quantity  in  parenthesis  is  measured  at  the  inlet  of 

aeration  basin 

o  =  quantity  in  parenthesis  is  measured  at  the  outlet  of 

aeration  basin 


C.   =  $100.00  for  1  <   P  <  20  horsepower 


C   =  $  21.57  for  20  <  P  <   200  horsepower 


1.00  $/sq.  ft. 


C3  =  9.259  x  10"3     $/cu.  ft. 
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4 

3.333  x  10"2 

s  - 

6.764 

C6  " 

5000. 

C7  ■ 

0.02 

C8  - 

1492. 

E    = 

0.75 

I 

k  = 

0.75 

M  « 

20 

Q  - 

45,120 

r  = 

0.10 

R  ■ 
o 

3.2 

T  - 

8760 

$/sq.  ft. 
S/sq.  ft. 

$ 

$/horsepower-hour 

$/horsepower 


(hour) 


years 

cu.    ft. /hour 

lb.  oxygen/horsepower 
hours 


x:  m      (BOD).  =  0.0109025  lb.  oxygen/cu.  ft.  solutic 


Determination  of  the  aerator  motor  size  of  the  inlet  section  of  the 
aeration  basin  (that  portion  of  the  lagoon  used  for  directly  oxidizing  the 
idealized  inorganic  component)  was  not  considered  to  be  part  of  the  opti- 
mization problem  since  the  size  is  fixed  once  the  waste  water  flow  rate 
and  the  BOD's  and  COD's  of  the  aeration  basin  inlet  and  outlet  streams 
are  given.   The  size  of  the  inlet  aerator  motor,  P  ,  was  calculated  using 
Equation  (27)  and  the  numerical  values  of  the  constants  given  previously. 


R  E 


l(C0D)i  -  (BOD).]  -   f  (C0D)q  -  (B0D)o  ]]    (27) 


A  value  of  220  horsepower  was  obtained  for  P  . 
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The  outlet  waste  solution  BOD  was  used  as  a  parameter  in  the  optimiza- 
tion study.   Three  cases  were  considered  in  which  it  was  required  to  obtain 
reductions  of  the  inlet  BOD  of  76,  90,  and  99  percent  respectively.   For 
the  76  percent  inlet  BOD  reduction  case,  optimum  policies  were  calculated 
for  one,  two,  and  three  stage  processes.   The  least  expensive  process  is 
the  one  stage  process.   Optimum  policies  were  calculated  for  processes 
composed  of  from  one  to  four  stages  for  the  90%  BOD  reduction  case.   For 
this  case,  the  two  stage  process  is  the  least  expensive  process.   Similarly, 
optimum  policies  for  processes  with  stages  numbering  from  one  to  six  were 
calculated  for  the  99  percent  BOD  reduction  case.   For  this  case,  the  four 
stage  process  is  least  expensive.   The  results  for  the  cases  of  76,  90,  and 
99  percent  inlet  BOD  reduction  are  summarized  in  Tables  1,  2,  and  3,  respective 
ly.   The  optimum  process  for  each  case  is  given  in  Table  4. 

Since  the  maximum  principle  does  not  guarantee  a  minimum,  it  is  necessary 
to  do  simulation  studies  with  a  mathematical  model  of  the  lagoon  to  insure 
that  a  minimum  has  been  found.   Simulation  studies  were  performed  with  com- 
puter program  SIM  which  was  programmed  to  simulate  the  process  taking  place 
in  the  lagoon  aeration  basin  based  upon  performance  Equations  (22)  through 
(25). 

The  results  of  the  simulation  studies  for  the  99  percent  BOD  reduction 
case  are  given  in  Table  5.   Three  cases  given  in  Table  5  are  slightly  less 
expensive  than  the  optimum  process  found  by  the  maximum  principle  solution. 
These  less  expensive  processes  resulted  because  of  discretization  error 
which  was  introduced  by  computer  program  OPT  during  the  optimization 
calculations.   The  optimum  process  can  be  computed  more  accurately,  however, 
the  improvement  obtained  in  the  total  process  cost  will  not  be  significant 
enough  to  justify  additional  computations.   In  the  fourth  simulation  study 
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given  in  Table  5,  the  aerator  motor  size  of  each  stage  was  obtained  by 
rounding  off  the  optimum  aerator  motor  sizes  of  the  first  three  stages 
to  integer  sizes  as  would  be  required  should  it  be  desired  to  purchase 
the  aerator  motors.   The  total  process  cost  of  the  fourth  simulation 
study  is  not  significantly  more  than  the  total  cost  of  the  optimum 
process . 

For  the  model  used  in  this  study,  several  trends  in  the  relationships 
among  the  variables  are  apparent  upon  inspection  of  Tables  1,  2,  3,  and  4. 
At  a  given  conversion,  the  optimum  total  volume  required  to  achieve  this 
conversion  decreases  as  the  number  of  stages  increases.   Both  the  total 
cost  and  the  optimum  number  of  stages  for  the  BOD  reduction  section  of  the 
lagoon  increases  as  the  percentage  BOD  conversion  increases.   For  all  of 
the  multi-stage  processes,  tapered  aeration  as  discussed  by  Sawyer  (8), 
is  the  best  policy.   This  is  due  to  the  fact  that  the  rate  of  BOD  reduc- 
tion decreases  as  the  treatment  progresses. 

The  discrete  maximum  principle  shows  promise  as  a  useful  technique 
for  the  optimization  of  aerated  lagoons  which  are  operated  in  a  stage- 
wise  fashion.   The  use  of  this  technique  with  the  aerated  lagoon  model 
used  in  this  work  predicts,  as  was  expected,  that  tapered  aeration  is  the 
optimal  operating  policy  for  multi-stage  lagoons. 

For  all  of  the  optimal  cases,  the  total  aeration  volume  in  each  case 
is  an  order  of  magnitude  smaller  than  the  Sugar  Creek  aeration  basin  volume 
(about  1,700,000  cu.  ft.).   This  may  be  partially  due  to  the  value  of  the 
first  order  reaction  rate  constant,  k  ,  which  was  used  in  the  numerical 
computations.   Too  large  a  value  of  k  would  account  for  the  smaller 
optimal  aeration  volume  calculated.   An  improved  model  with  calculations 
to  support  this  belief  is  given  below. 
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An  improved  lagoon  model  can  be  derived  by  making  material  balances 
for  both  the  microorganisms  and  the  BOD  component  around  each  of  the 
ideal  backmix  stages  of  the  lagoon.   The  material  balances  for  the  micro- 
organisms and  the  BOD  component  are  given  by  Equations  (28)  and  (29), 
respectively. 


„n-l 


n-1 


C  '  + 


kM  xn  rn 
"  Xl  CB 

K  +  x? 


kM  x11  Cn 
R  Xl  S 

Y(K  +  x") 


h 


(28). 


(29) 


where 


L       ■  the  microorganism  concentration  in  the  nth  stage 
k       =     the  endogeneous  respiration  rate  constant 
8n  =  Vn/Q 
If  the  power  requirement  for  each  stage  is  calculated  by  Equation  (9) 
and  if  the  total  stage  cost  is  calculated  by  Equation  (16),  then  for  a  two 
stage  process  a  one  dimensional  search  can  be  used  to  determine  the  values 
of  P   and  V   which  will  give  the  lowest  aeration  basin  cost. 

A  one  dimensional  exhaustive  search  was  performed  for  76,  90  and  99 

percent  BOD  reduction.   The  search  variable  used  was  x   ,  the  first  stage 

1  1    ?   1 

outlet  BOD.   With  x.   known,  and  the  constants  given  below,  C_  ,  C_  ,  8  , 
1  B    B 

2    12 
8  ,  P  ,  P  and  the  total  cost  of  the  process  may  be  calculated.   The  ex- 
haustive search  procedure  was  used  to  calculate  the  total  process  cost  at 

1  0       2 

a  series  of  x^  values  taken  at  intervals  between  x.      and  x.  .   The  process 

with  the  xl     value  which  gave  the  lowest  total  cost  was  selected  as  the 

optimum.   The  above  procedure  was  programmed  for  use  on  an  IBM  1410  computer. 
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The  name  of  the  computer  program  which  was  used  to  perform  the  one  dimen- 
sional search  calculations  is  SEARCH.   Program  SEARCH  is  described  in 
Appendix  III.   The  numerical  values  of  the  constants  used  in  program 
SEARCH  are  as  follows: 


Cfi  «  0.0  ppm  K  =  100  ppm 

kM  =  0.1  hr"1  x°  =  175  ppm 

y.  m  0.002  hr"*1  Y  =  0  5  lt)s  ■  microorganisms  formed 

«  *       lb.  BOD  consumed 


The  optimum  results  obtained  from  the  exhaustive  search  for  the  three 
cases  are  given  in  Table  6. 

For  the  case  of  76  per  cent  BOD  reduction  a  one  stage  process  was 
optimum.   If  the  values  of  C   ,  k,  K,  and  x.   for  the  optimum  76  per  cent 
BOD  reduction  case  are  used  in  the  equation  which  defines  k',  a  value  of 
k  equal  to  0.0437  results.  This  new  k  value  when  used  in  computer 
program  SIM,  the  lagoon  process  simulation  program,  gives  a  lagoon  volume 
of  3,272,449.  cu.  ft.  and  a  total  cost  of  $421,824.  for  the  one  stage  76 
per  cent  BOD  reduction  case.   This  lagoon  volume  is  of  the  same  order  of 
magnitude  as  the  one  at  Sugar  Creek. 

The  minimum  cost  solution  for  90  per  cent  BOD  reduction  is  a  two  stage 
process.   Using  this  solution  and  the  defining  equation  for  k  ,  values  of 
0.0433  and  0.0619  were  obtained  for  the  first  and  second  stages,  respectively. 
These  numbers  indicate  that  the  assumption  of  a  constant  k'  value  throughout 
the  process  is  not  reasonable.   In  the  next  section  of  this  report,  a  lagoon 
aeriation  basin  model  is  proposed  which  gives  a  more  detailed  account  of 
the  processes  which  occur  in  the  aeration  basin.   Optimization  of  this 
model  should  give  more  realistic  answers. 
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5.0   PROPOSED  FUTURE  WORK 

More  accurate  performance  equations  for  the  aerated  lagoon  process 
can  be  obtained  by  quantitatively  considering  the  primary  components 
which  are  present.   The  model  used  in  the  present  optimization  study 
assumed  that  the  BOD  reduction  rate  could  be  approximated  by  a  pseudo 
first  order  kinetic  expression  written  in  terms  of  BOD  concentration.   As 
shown  in  the  last  section  this  assumption  is  not  valid  since  the  pseudo 
rate  constant,  k  ,  varies  appreciably  throughout  the  process.   Further- 
more, the  effects  of  oxygen  mass  transfer  rate  and  oxygen  concentration 
upon  reaction  rates  need  to  be  considered  in  detail. 

The  kinetic  and  mass  transfer  expressions  for  the  aerated  lagoon 
model  which  is  proposed  for  future  work  are  derived  in  the  following 
paragraphs. 

Monod  (9)  first  used  a  model  to  express  the  rate  of  growth  of  micro- 
organisms which  represents  the  rate  as  being  proportional  to  the  nutrient 
or  substrate  concentration  at  low  nutrient  concentrations  but  reaching  a 
limiting  value  at  high  nutrient  concentrations: 
n 

xi 

K  +  x. 


where, 


r       =     rate   of  microorganism  growth 


M 
k       «     growth  rate  constant 


c.      =     nutrient   concentration 


«„     =     microorganism  concentration 


Michaelis-Menten  constant 
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This  model  assumes  that  sufficient  oxygen  is  present  to  satisfy  the 
needs  of  the  organisms. 

The  problem  of  maintaining  adequate  oxygen  concentrations  has  been 
recognized.   However,  little  work  has  been  done  to  determine  adequate 
oxygen  levels  even  for  pure  cultures. 

As  a  means  of  taking  into  account  the  effect  of  oxygen  concentration 

upon  the  rate  of  formation  of  organisms,  a  modification  of  Equation  (30) 

is  used  in  this  work. 

n  n 

Xl  X3 

r   =  kMx^   ( — )  ( — )  (31) 

B         2   K  +  x°  KQ  +  x^ 

where, 

x   ■  oxygen  concentration 

K_  ■  a  constant  similar  to  the  Michaelis-Menten  constant 

In  aeration  systems,  microorganisms  are  lost  through  a  process  known 
as  endogenous  respiration.  In  effect,  they  are  consumed  or  oxidized  by 
other  microorganisms.  McKinney  (10)  used  a  first  order  kinetic  equation 
to  express  the  endogenous  respiration  rate,  r  .   The  equation  can  be 
written  as 

re  "  hX2  (32) 

where  k_  is  the  endogenous  respiration  rate  constant. 

By  analogy  to  Equation  (31),  the  factor,  (x„/(KQ  +  x  ))  , 
is  used  in  this  work  to  take  into  account  the  effect  of  low  oxygen 
concentrations  upon  the  endogenous  respiration  rate 


\>   *2  ( " >  <33) 


h+% 
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If  it  is  assumed  that  the  microorganism  population  increase  is  pro- 
portional to  the  decrease  in  nutrient  concentration,  then  the  nutrient 
reduction  rate,  r   ,  can  be  given  by 

rx  =  "  —  C34) 

where  Y  is  the  yield  factor. 

The  oxygen  consumption  rate,  r„,  is  equal  to  the  sum  of  the  endogenous 
respiration  and  nutrient  reduction  rates 

rn  =  r  +  r  (35) 

0      ex  K      J 

Oxygen  is  supplied  to  a  stage  of  the  aeration  basin  by  the  aerator. 

Since  aeration  is  a  mass  transfer  operation,  it  is  taken  into  account  by 

letting  the  mass  transfer  rate  per  unit  of  stage  volume  be  equal  to  the 

product  of  a  combined  overall  liquid  phase  mass  transfer  coefficient  and 

interfacial  area,  ICa,  and  a  driving  force  which  is  given  by  the  difference 

between  the  oxygen  saturation  concentration  and  the  oxygen  concentration  in 

the  stage,  (x,  -  x„).   To  calculate  oxygen  input  rate  as  a  function  of 

aerator  horsepower,  the  mass  transfer  coefficient  is  written  as  a  function 

of  the  aerator  horsepower,  shaft  speed,  and  the  stage  volume  as  follows: 

Kja  =  Kc  (£l-)C  (vVb  (36) 


where 


and 


b,  c,  and  K  are  constants 


P  =  aerator  horsepower 

R  ■  aerator  shaft  speed 

Trn  _ 

V  =  stage  volume 
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Equation   (36)  was   suggested  by   the   investigations   of  Blakebrough 
and   Sambamurthy   (11),   which  were   conducted   in  batch  laboratory   fermentors. 

A  nutrient,    or  BOD,   mass   balance   for  an  ideal  backmix  stage  of  the 
nutrient  degradation  section  of  the   lagoon  may  be  written  as 

j   n 

OX,  , 

„n        1  .,     n-l    ,    TTn  ,     n  /„~,^ 

V     =     Q  x         +  V     r     -   Q  x  (37) 

dt  L  xi 

The  terms  of  Equation  (37)  take  into  account  accumulation  of  nutrient 
in  the  stage,  nutrient  washing  into  the  stage,  nutrient  degradation  in 
the  stage,  and  nutrient  washing  out  of  the  stage,  respectively. 

Similar  mass  balances  for  microorganisms  and  oxygen  are  given  by 
Equations  (38)  and  (39),  respectively. 


dx„  , 

V  =  Q  x    +  V  r  -  Q  x  (38) 

dt        2         B      2 


dx,         , 
V^  -  Qx^  +  Vn 


Kj.a  (x3  -  x3)  +  rQ 


Q  *3     (39) 


The  performance  equations  for  the  model  which  is  proposed  for  future 
work  are  summarized  in  Table  7.   When  used  in  conjunction  with  a  stage 
cost  equation  such  as  Equation  (23),  this  model  can  be  used  to  obtain 
more  accurate  values  of  the  independent  variables  of  the  optimization 
study. 
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7.0  NOMENCLATURE 

al  al 

Constants  used  in  Equation  (10) 

b 

Constant  used  in  Equation  (30) 

c 

Constant  used  in  Equation  (30) 

BOD 

Biological  oxygen  demand 

S 

Aerobic  microorganism  concentration 

«$l 

Initial  investment  cost  of  aerator  motor 

(C»)2 

Cost  of  aeration  basin 

(Cg) 

Present  worth  of  electrical  power  cost  which 
is  paid  over  a  period  of  years 

COD 

Chemical  oxygen  demand 

Cl>  Ci 

Constants  used  in  Equation  (10) 

C2 

Unit  land  cost  of  aeration  basin 

S 

Unit  digging  cost  of  aeration  basin 

C4 
S 

Unit  cost  for  rock  lining 

2 
Constant  defined  as  — 

of  aeration  basin 
r       2]  2/3 
ff- J        (6C2  +  7C4) 

C6 

Initial  cost  of  aerator  turbine  and  supporting  stand 

C7 

Electrical  power  cost 

C8 

Constant  defined  as 

M 
C7  T  Z      (1  +  r)"1 
1-1 

E 

Mechanical  efficiency  of  aeration  unit 

f  (x,    ;  S   )   Right  hand  side  of  optimal  recurrence  relation, 

Equation  (26) 

142 


g  (x^;  6   )  Left  hand  side  of  optimal  recurrence  reaction, 

Equation  (26) 

G  ,G(x,   ;  9  )  Total  cost  of  n-th  aeration  basin 

h  Aeration  basin  depth 

k  Pseudo  first  order  reaction  rate  constant 

K  Microorganism  growth  rate  constant 

k-  Endogenous  respiration  rate  constant 

K  Michaelis  -  Menton  constant 

Kc  Constant  used  in  Equation  (30) 

Kra  Combined  overall  liquid  phase  mass  transfer 

and  interfacial  area 
Kq  Constant  similar  to  K 

-<  Length  of  aeration  basin 

M  Useful  life  of  aeration  basin 

n  Superscript  denoting  the  n-th  stage 

of  the  aeration  basin 
N  Superscript  denoting  the  last  stage  of  the 

aeration  basin 
P  n-th  stage  aerator  motor  size 

P  Aerator  motor  size*  for  inlet  section  of  aeration 

basin  which  directly  oxidizes  the  idealized 
inorganic  component. 
Q  Waste  water  flow  rate  to  aeration  basin 

r  Simple  interest  rate 

£B  Aerobic  microorganism  formation  rate 

r  Endogenous  respiration  rate 
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rQ  Oxygen  utilization  rate 

r  BOD,  i.e.,  nutrient  reduction  rate 

R  Aerator  shaft  speed 

R  Oxygen  transfer  rate  constant 

S  Total  cost  of  BOD  degradation  section  at  aeration  basin 

t  Time 

T  Aeration  basin  operating  hours  per  year 

V  Volume  of  n-th  aeration  basin  stage 
w  Width  of  lagoon 

x  BOD  concentration 

x   ,  x  BOD  concentration  in  (n-l)-th  and  n-th  stage  of 

aeration  basin,  respectively 

N 
x.  BOD  concentration  in  waste  water  leaving 

aeration  basin 

0 
x,  BOD  concentration  in  waste  water  flowing  into 

aeration  basin 
*2  Sum  °f  costs  of  all  stages  of  the  aeration  basin 

up  to  and  including  the  n-th  stage. 

Used  in  Equations  (21)  and  (24) 
*2  Microorganism  concentration  in  n-th  stage  of 

aeration  basin.   Used  in  Equation  (31)  and  thereafter. 

*2  Oxygen  concentration  in  n-th  stage  of  aeration  basin. 

s 

x.  Oxygen  saturation  concentration 

Y  Yield  constant  /  lbs,  microorganisms  formed* 

lb.   BOD  consumed         ' 
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Greek  Letters 

Oi. 

The  power  to  which  a  term  in  Equation  (23)  is  raised 

f>  The  coefficient  of  a  term  used  in  Equation  (23) 

9  Volume  of  n  th  stage  as  used  in  Equations  (22)  through  (24) 

9  Holding  time  of  n  th  stage  of  aeration  basin  (Vn/Q).  Used 

in  Equations  (28)  and  (29). 
9  Aerator  motor  size.   Used  in  Appendix  I. 
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APPENDIX  I 
EXPLANATION  OF  LAGOON  OPTIMIZATION  COMPUTER  PROGRAM  - 
OPT  AND  SUBROUTINE  -  DMP 

To  simplify  programming  of  OPT,  the  computer  program  which  was  used 
to  perform  the  optimization  calculations,  different  performance  equations 
were  used  in  place  of  the  performance  equations  given  in  Section  3.2.   The 
transformation  equation  which  was  used  in  OPT  is  obtained  by  rearranging 
Equation  (9)  and  defining  the  electric  motor  size  as  the  n  th  stage 
decision  variable,  9  . 


x*-1   -  9n  RoE/Q  =  TCx*'1;  Sn)  (W) 


where 


xl     ~     oxygen  demand  of  the  waste  solution  flowing  into 

stage  n 
xl     =  oxygen  demand  of  the  waste  water  solution  flowing  out 

of  stage  n 


e" 


aerator  electric  motor  size  at  stage  n 


RQ  =  oxygen  transfer  rate  constant 

E   =  mechanical  efficiency  of  aeration  unit 

Q   ■  waste  water  flow  rate  into  aeration  basin 

The  stage  volume  required  for  stage  n,  Vn,  is  calculated  by  combining 
Equations  (8)  and  (9)  so  that  the  outlet  BOD  is  eliminated  from  the 
resulting  equation  which  is  given  below. 

_2_  " 


»S1     -  6n} 

Re 


The  cost  of  stage  n  is  calculated  using  Equation  (1-3)  which 
is  obtained  by  substituting  Equation  (1-2)  into  Equation  (16). 


where 


Bfr£"l|  9n) 


C6  +  C8  9n  +  /?(9n)   +  C3 


+  C. 


R  E 


°<  =  a  when  1  <  9n  <  20 


0 


Cx   when  1  £  8  <  20 


ft    <•     C^  when  20  <  9n  <   200 
The  accumulated  cost  of  the  process  is 


0  9' 


i  0*1 
k  (— — 


R  E 
o 


(1-3) 
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n      n-1     .  n-1  _n,    0 
<2     ■  x2   +  G(x   ;  9  )  ,  x 


(1-4) 


The  objective  function  which  is  to  be  minimized  is 


S  =  x ' 


(1-5) 


The  transformation  equation,  Equation  (1-1),  may  be  used  to  calculate 
the  BOD  concentration  at  the  outlet  of  the  nth  stage  given  the  inlet  BOD 
concentration  and  the  electric  motor  size  of  the  stage  (the  decision 
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variable).   The  volume  of  stage  n  may  be  calculated  using  Equation  (1-2). 
The  cost  of  stage  n  may  be  calculated  using  Equation  (1-3) .   Equation  (1-4) 
gives  the  total  cost  of  all  stages  up  to  and  including  stage  n.   The  total 
cost  of  the  BOD  reduction  section  of  the  lagoon  is  given  by  Equation  (1-5). 
Equations  (1-1)  through  (1-5)  are  the  performance  equations  which  were  used 
in  programming  computed  program  OPT. 

For  the  transformation  and  stage  cost  equations  given  in  this  Appendix, 
the  left  hand  side  of  the  necessary  optimal  recurrence  equation,  Equation 
(26) ,  is  given  by  Equation  (i-6) .   The  right  hand  side  of  Equation  (26)  is 
given  by  Equation  (1-7) . 


,  n-1   „n. 

;(x1  ;  e  ) 


S  +  a2C5 


JL-JL 


Qx 


n-1 


k'(- 


,,  n  „n+l. 
fCx^  6    ) 


C3  +  a2C5 


Q  1 


n+l 


a2-l 


k'(S-e-) 


„-I 


+  aB(en+1)a_1  +  c0 


+  dUe")0-1  +  c0 


(1-6) 


.n+l 


(1 


(1-7) 


where 


for   1  <  8n  <  20   or   1  <  en+1  <  20 


for  20  <  8n  <  200  or  20  <  en+1  <  200 


for   1  < 


<  20   or   1  <  8n+1  <  20 


=  C, 


for  20  <  Qn  <  200   or  20  <  6n+1  <  200 
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The  values  of  o(.  and  /3   used  in  this  problem  change  abruptly  at 

~n      ~n+l 

W   or  9    equal  to  twenty  because  the  aerator  motor  cost  versus 

size  relation,  which  was  used  in  this  optimization  study,  has  a  dis- 
continuity at  a  value  of  twenty  horsepower.   This  abrupt  change  in 

CL   and  /3    causes  the  functions  g  (x°  }  9n)  and  f  (x";  9n+1),  when 
plotted  versus  9  and  9   ,  respectively,  to  be  discontinuous  at  9n  or 
9    equal  to  twenty.   The  discontinuity  in  the  f  (x"j  9   )  versus  9n+1 
curve  is  shown  in  Figure  1-1  for  a  typical  f  (x.n;  9   )  function.   Even 
though  the  derivation  of  the  necessary  recurrence  equation,  Equation  (26), 
requires  that  the  functions  g  (s^   |  9n)  and  f  (x  n;  8n+1)  be  continuous 
over  the  range  of  8  and  8    values  encountered  in  the  process,  it  was 
found  that  the  maximum  principle  may  still  be  used  when  the  objective 
function  is  piecewise  smooth.   The  necessary  recurrence  equation  is  valid 
over  each  piecewise  smooth  segment  of  the  f  (x.n;  9   )  versus  8nTl  curve. 
The  search  technique  which  is  used  to  find  9n   is  modified  to  take  into 
account  the  discontinuity.   This  procedure  allows  the  search  technique  to 
search  the  whole  range  of  9    values.   The  search  technique  and  the  modi- 
fication which  was  used  is  discussed  in  the  following  paragraphs  which 
describe  the  computer  program  OPT. 

Computer  program  OPT  can  be  easily  described  by  reviewing  the  compu- 
tational procedure  given  in  Section  3.3.   By  starting  at  the  first  stage 
of  the  process  and  assuming  a  value  of  9  ,  a  fixed  value  of  g  (x  °;  81) 
may  be  calculated  using  Equation  (1-6).   The  BOD  at  the  outlet  of  stage 
one,  x^  ,  may  be  calculated  by  using  the  transformation  equation,  Equation 
(1-1).   A  value  of  8  may  be  calculated  by  using  the  necessary  recurrence 
equation  which,  for  this  problem,  is  obtained  by  equating  Equations  (1-6) 
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and  (1-7).   The  necessary  recurrence  equation  has  the  same  form  as 
Equation  (26a). 

Due  to  the  complexity  of  the  right  hand  side  of  the  necessary  re- 
currence equation,  i.e.,  the  function  f  (x  n;  9n   ),  it  is  expedient  to 

2 
use  a  search  technique  to  solve  for  6  .   Similarly,  a  value  of  the  de- 
cision variable  for  each  stage  of  the  process  may  be  calculated  by 
alternate,  repeated  applications  of  the  transformation  equation  and 

searching  the  necessary  recurrence  equation.   If  the  value  of  the  outlet 

N 
BOD,  x.,  is  not  equal  to  the  required  value,  the  above  procedure  is 

repeated  by  assuming  a  new  value  for  9  .   If  the  value  of  x,  is  suffici- 
ently close  to  the  required  value,  no  further  calculations  are  required. 

Computer  program  OPT  was  programmed  to  carry  out  the  computational 
procedure  given  above.   When  work  was  first  commensed  on  OPT,  it  was 
assumed  that  only  one  value  of  9    would  be  found  when  searching  the 
necessary  recurrence  equation  over  the  range  of  9°   values  which  are 
feasible  for  the  process.   This  situation  is  indicated  by  the  dashed  line 
in  Figure  1-1.   Several  trials  using  the  Bolzano  search  technique,  which 
is  described  by  Stanton  (1),  revealed  that  the  computation  time  for  this 
technique  was  not  excessive.   Consequently,  the  Bolzano  technique  was 
selected  for  use  in  OPT  to  search  for  9n 

After  the  Bolzano  technique  was  incorporated  into  OPT,  very  small 
values  of  the  state  variables  were  obtained  for  some  cases  tried.  Further 
investigation  revealed  that  two  values  of  9n   satisfy  the  necessary  re- 
currence equation.   This  situation  is  shown  by  the  solid  curve  for 

,,  n  „n+l. 

IVX^;  w   ;  on  iigure  1-1.   Optimum  processes  were  not  obtained  using  the 

smaller  9    values.   Therefore,  the  Bolzano  technique  was  made  to  converge 

to  the  larger  9    values.   Convergence  to  the  larger  9n   values  was 

assured  by  starting  the  Bolzano  technique  on  the  portion  of  the  f(xn-  9n+1) 
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versus  9    curve  which  has  a  positive  slope. 

It  was  discovered  that  when  the  9    value  being  searched  for  was 
near  the  discontinuous  point  of  the  f(x";  9n  )  versus  9n   curve,  the 
Bolzano  technique  failed  to  converge  to  a  9    value.   The  reason  for 
the  failure  of  the  Bolzano  technique  to  converge  to  a  9    value  was 
not  investigated  because  a  simple  solution  for  this  difficulty  was  found. 
It  was  observed  that  the  Bolzano  technique  converged  to  a  8n   value  at 
positions  away  from  the  discontinuity  in  only  a  few  trials.   Consequently, 
a  limit  was  placed  on  the  number  of  trials  which  were  permitted  for  the 
Bolzano  technique  to  search  for  a  8   .In  the  event  that  the  number  of 
trials  used  by  the  Bolzano  technique  exceeded  the  limit,  i.e.,  in  the 
event  that  the  technique  failed  to  converge,  9  '  was  set  equal  to  the 
8    value  just  to  the  left  hand  side  of  the  discontinuity  because  this 
8    value  gives  a  lower  aerator  motor  cost  than  the  9n   value  just  to 
the  right  of  the  discontinuity. 

This  work  has  shown  that  it  is  possible  to  extend  the  use  of  the 
discrete  maximum  principle  to  solve  optimization  problems  which  have 
piecewise  continuous  objective  functions.  Rosenbrock  and  Storey  (2)  have 
pointed  out  that  piecewise  continuous  cost  versus  size  curves  are  frequently 
encountered  when  dealing  with  power  equipment. 

Subroutine  DMP  is  used  with  program  OPT.   Subroutine  DMP  consists  of 
several  defined  functions  and  it  is  used  to  reduce  the  size  of  OPT. 

Table  1-1  explains  the  symbols  used  in  OPT  and  DMP.  Logic  diagrams 
for  OPT  and  DMP  are  given  in  Figures  1-2  and  1-3,  respectively.  Program 
OPT  and  subroutine  DMP  are  listed  in  Table  1-2. 
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Table  1-1.   Explanat: 
for  Program 

Lon  of  Computer  Program  Variables 
OPT  and  Subroutine  DMP 

Symbol 

Explanation 

A 

Current  sum  of  G(x. ;  9   )  values, 

Al 

a.  on  a.' 

A2 

a2 

AC1 

C6 

AC2 

C8 

AC3 

Cl 

AC4 

C5 

AC5 

s 

B 

Current  sum  of  P  values 

Bl 

Q/RoE 

B2 

Q/k' 

C 

Current  sum  of  V  values 

DGDT 

Jc(x^;  9n+1)  /J  9n+1 

DGDX 

c)G(x";  9n+1)  /dxj 

DT 

.  Incremental  change  in  9 

DTDT 

<)T(x^  9n+1)  /  cl9n+1 

DTDX 

c)T(x*;  9n+1)  /  c)x^ 

DUDX 

c)Vn  /  o>  x^ 

D2GDT 

<)2G(x*;  9n+1)  /  (d9n+1)2 

D2UDT 

c)V  /  (t)  9n+1)2 

E 

Allowable  error 

EFF 

E 

F 

A  defined  function  used  to  simplify  the 
computer  program 
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Table  1-1.  (con't) 


Symbol  Explanation 

G  G<X*J  en+1) 

GL  Left  hand  side  of  optimal  recurrence 

equation 

GR  Right  hand  side  of  optimal  recurrence 

equation 

G1,G2,G3         Defined  functions  used  to  simplify 
calculation  at  D2GDT 

I  A  counter  used  as  the  stage  number 

II  A  counter  used  to  determine  the  number 

at  times  the  Bolzano  technique  is  applied 
at  each  stage. 

IT  Maximum  number  at  times  Bolzano  technique 

is  applied  at  each  stage 

J  A  counter  used  in  the  logic  steps  of  the 

Bolzano  technique 

K  Number  of  stages  in  process 

Q  Q 

QA  Multiplication  factor  used  to  obtain  a 

positive  value  of  D2GDT. 

RK  k' 

S  A  dummy  variable  used  to  sum  G(x1 ;  8   ) 

ST  A  dummy  variable  used  to  sum  P 

SU  A  dummy  variable  used  to  sum  V 

T  0n+1 

TLIM  Maximum  value  of  9    that  can  be  used  in 

remaining  stages 

n+1 
Tl  Initial  value  of  9    for  stage  1. 


U  Stage  volume 


Table  1-1.   (con't) 


X 
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Symbol  Explanation 


WT  R. 

o 


XA  A  dummy  variable  used  to  calculate  the 

stage  outlet  BOD  concentration  in 
program  SIM 

XF  x° 

XO  Dummy  variable  used  with  x. 

(1 


N 
XR  Desired  X,  value 


X1,X2,X3  Values  of  8    during  Bolzano  search 

Y1,Y2,Y3  Values  of  right  hand  side  of  optimal 

recursion  equation  during  Bolzano  search 

Z  Value  of  left  hand  side  of  optimal 

recursion  equation  during  Bolzano 
search 


155 


«       J       3 


156 


ENTER   SUBROUTINE 
DMPCX.T,  U,G,6L,GR,D2G0T,  A2,Bl,  B2,  ACI,  AC2.AC4,  AC5) 


NO>-»- 


-e-frEs) 


AC3»2I.57 
Al'1.08 


AC3«  100. 
Al  «  0.53 


RETURN 


Fig.  1-3..   Logic  diagram  for  subroutine  DMP. 


KONSS      JC8 

MSNSS      COMT  10  MINUTES.  30  PAGES. 

MONSJ      ASGN  MJB.12 

MCN$$      ASGN  MG0.16 

MCNSS      MODE  GO 

MCNSI  EXE«    FORTRAN. ..... .CPT 

1  FORMAT (7E 13. 6) 

2  FORMAT! I  3 . 3X , F9 . 6 , 3X , F6 . 1 . 3X . F9.0 . 3X . F9. 0 ) 

3  FORMAT(2E13.6,2I3) 

4  F0RMAT(6X,F9.6) 

5  FORMAT (13X,F6.1,3X,F9.0,3X,F9.0) 
RFADd.llQ.XF.XR 

READ! 1,1) AC1.AC2.AC4.AC5 

READ11 ,1) EFF.RK.WT.A2 

READ (1.1) E 

WRITE! 3.1 10, XF.XR 

WRITEI3.1 )AC1,AC2,AC4,AC5 

WRITE(3,1 1EFF.RK.WT.A2 

WRITE(3.1)E 
14    RFADf 1,3) Tl.DT.K.IT 

WRITF(3,3)T1.DT,<, IT 

XO  =  XF 

1=1 

J  =  l 

I  1  =  0 

31=Q/(WT*EFF) 

32=0/RK 
22    T=T1 

X  =  XG 

CALLnMP(X,T.U,G,GL.GR,D2GDT.A2.31.32.ACl,AC2,AC4.AC5) 
WRITEI3.2) I .X.T.U.G 

2  6    CA!_LDMP(X,T,U,G,GL.GR.D2GDT.A2.B1,82.AC1.AC2.AC4.AC5> 

Z  =  GL 

!F(  I.NE.DGOTC35 
29    S=A+G 

A  =  S 

ST  =  5-i-T 

8=ST 

SU=C+U 

C  =  SU 
35  X=X-T/31 

TLIM=81*X 

0A=.5 
38  T=QA*TLIM 

3  9  CALLDMP(X,T,U,G,GL.GR.D2GDT,A2.B1.32»AC1.AC2.AC4.AC5) 

IFU.NE.l  JG0T053 
IF(D2GDT.LT.O. 1G0T046 
IFIABSi (2-GR)/Z).LT.E)GOT070 
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'.Cable    1-2.     (Con't) 

IFI6R.LT. Z)S0TC48 
44   QA=OA-.l 

GCTC38 
46    QA=QA+.05 

GCTC38 
48  Y1=GR 

X1  =  T 

J-2 

T=.09*TLIM 

GCTC39 
53  Y2  =  GR 

X2  =  T 
55  II  =  H  +  1 

T=.5*(X1+X2> 

CALLDMP(X,T»U,G,GL.GR,D2GDT.A2,B1,B2,AC1,AC2,AC4,AC5) 

Y3  =  GR 

X3  =  T 

IF(  I  I.GT. IT1GCTC69 
IF!ABS( (Y3-Z)/Z).LT.E)GCTC70 
IF(Y3.GE.Z)GCT»66 
63  Y1  =  Y3 

X1  =  X3 

GCT055 
66  Y2=Y^ 

X2  =  X3 

GCTC55 
69  T=19.9 
7  0  1=1+1 

S=A+G 

A  =  S 

ST=3+T 

B=ST 

SU=C+U 

C  =  SU 

WRITE13.2) I ,X»T,U»G 

11=0 

IF(  I.LT.OGCT026 
81  X=X-T/31 

WRITE(3.4)X 

WRITE13.5 1B.C.A 

1=1 

J  =  l 

A=.0 

B  =  .0 

C  =  .0 

IF(X.LT.XR)GCTC14 
90  T1=T1+DT 

GCTG22 

END 
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Table    1-2.    (Con't) 


MCNSS  EXEQ    FORTRAN* »»»»» . 

SUBROUTINE    DMP ( X .T .U .G .GL.GR.D2GDT » A2 .31 »B2 »AC1 .AC2 . AC4 . AC5 ) 
1    FCRMAT(7E13.6) 

IF(T.LE.20.)GCTC6 
3  AC3»21-57 

A1=1.0B 

GCT08 
6    AC3=100. 

Al=0.53 
8    F«81*82/« (B1*X-T)*(B1*X-T) ) 

DTDT=F»X 

U=B2*T/(B1*X-T) 

G=AC1+AC2*T+AC3*T**A1+AC4*U**A2+AC5*U 

DG0T=AC2+A1*AC3*T**(A1-1.)+(AC5+A2*AC4*U**(A2-1.) )*DTDT 

DTDX=B1 

DUDX=-F*T 

DGDX=DUDX*( AC5*A2*AC4*U**( A2-1. ) I 

GL=DGDT 

GR=DGDT+DGDX/B1 

D2UDT=F*X*2./(B1*X-T) 

G1=(A1*A1-A1 )*AC3*(T**(Al-2. ) ) 

G2= (AC5+A2*AC4*U**( A2-1. ) )*D2UDT 

G3=< (A2*A2-A2)*AC4*(U)**(A2-2.) )*DTDT 

D2GDT=Gl+G2+G3 

RETURN 

END 
MCNSS      EXEQ  LINKLCAD 

CALL  CPT 
MCNSS      EXEQ  CPT.MJE 


DATA 

.451170E+05  .109025E-01  .109025E-02 

.500000E+04  .149158E+04  .676441E+01   .009260E+00 

.750000E+00  .750000E+00  .320000E+01   .666667E+00 
,0000505+00 

.1845C0E+03  .10CCOOE+00  2  20 

.139600E+03  .100000E+00  2  20 

.1C7700E+03  .100000E+00  3  20 

•  081200E+03  .lOOC'OOE+00  4  20 
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APPENDIX  II 
EXPLANATION  Of  LA30QN  SIMULATION  COMPUTES  paOGBAM  -  SIM 

Computer  program  SIM  simulates  the  operation  of  the  aer- 
ated lagoon  according  to  the  aerated  lagoon  model  used  In  the 
discrete  maximum  principle  optimization  study.  The  symbols 
used  in  program  SIM  are  the  same  as  those  given  in  Table  1-1. 
A  logic  diagram  for  SIM  is  given  in  Figure  II-l.  program  SIM 
is  listed  in  Table  II-l. 
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READ: 
Q.XF.XR 

ACI,AC2,AC4,AC5 
EFF,RK,WT,A2 


WRITE: 
Q,  XF,  XR 

ACI,AC2,AC4,AC5 
EFF,  RK,  V/T,  A2 


AC3  =  I00. 
Al  «0.53 


Fig.  II-l.   Logic  diagram  for  computer  program  SIM. 


MCNSS      JCS 

MCNSS    •   CCMT   5  MINUTESt  10  PAGES. 
KCNSS      ASG^  MJB.12 
"CNSS      ASG-  MG0,16 
MCNSS      N.COL  uC 

MCNSS      EXEQ  FORTRAN ,..,,, ,S IM 
1  FCRMATI7E13.6) 
.2  FORMAT! I  3 , 3X , F9. 6 , 3X ,F6 . 1 » 3X , F9 .0 , 3X .F9. 0 ) 

3  FORMAT; 13X.E13.6.I3) 

4  FORMAT (6X.F9. 6) 

5  FORMAT(I3X,F6.1,3X.F9.0,3X,F9.0) 
READ* lil )G,XF,XR 

P-ADd  >DACl»AC2iAC4iAC5 
READUil) EFF.RK.WT.A2 
WRITE! 3il)Q,XF»XR 
WRITE (3,1 ) AC1.AC2.AC4.AC5 
V.'R:TE(3,1)EFF,RK,WT.A2 
31=0/ ( WT*EFF ) 
32=CVRK 
14   READ(l>3JX»fC 

16    RFADil ,3) T 

■:F(T.LE.20.)GOT021 
18    AT3=21  .57 

Al-1.08 

G0T073 
21    AC3=100. 

Al=0.53 
23    U=32*T/(B1*X-T) 

G  =  AC:-rAC2*T  +  AC3*T**Al+AC4*U**A2+AC5*U 

S=A+G 

A=S 

ST  =  i'2  +  T 

B  =  ST 

SU=G+U 

C  =  SU 

WRITE (3, 2) I .X.T.U.G 

IF(  I.GE..OG0T036 

1  =  1+] 

X=X-T/S1 

G0T016 
36    WRITE<3»5 )3,C»A 

X=X-T/B1 

WRITE! 3.4 )X 

A  =  0. 

R  =  .C 

c-.o 

GC7014 
END 
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fable  II-l.  (Con't) 


MO,\SS 


MCNSS 


EXEQ  LINKLCAD 
CALL  SIM 
EXEO  SIM.MJB 


DATA 
.451170;+05 
.5C00C0E-r04 
.  750  0COE-i-00 


.109025E-01 
•14915SE+04 
.043700E-rOC 
•10902SE-01 
.155800E+03 


•261660E-02 
.676441E+01 
.320000E+01 

1 


.C09260E+00 
.666667E+00 
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APPENDIX  III 
EXPLANATION  OF  C0MPUT2H  P80GaAM  -  SEARCH 

Computer  program  SEA3CH  is  used  to  perform  a  one  dimen- 
sional exhaustive  search  on  the  two  stage  aerated  lagoon  model 
given  by  Equations  (28)  and  (29).   Table  III-l  explains  the 
symbols  used  in  SfiAtiCH  which  have  not  previously  been  explained 
in  Table  1-1.  The  logic  .diagram  for  SiiiuiCH  is  given  in 
Figure  III-l.   Program  SiitiriCH  is  listed  in  TaDle  III-2. 
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Table  III  -  1.   Explanation  of  the  Symbols  Used  In  Computer 
Program  SEARCH 


AK 
B 


Symbol  Explanation 


Constant  used  in  the  calculation  of  the 
incremental  change  in  the  search  variable 
X(2) 


kM 


Constant  used  in  the  calculation  of  the 
incremental  chanse  in  the  search  variable 
X(2) 

CK  K 

C(1),C(2),C(3)        Concentration  of  micro-organisms  going  into 
the  first  stage,  in  the  first  ideal  back- 
mix  stage,  and  in  the  second  ideal  back- 
mix  stage,  respectively. 

D  Constant  used  in  the  calculation  of  the 

incremental  change  in  the  search  variable 
X(2) 
DEL  Incremental  change  in  the  second  variable 

X(2) 

DK  H 

G(2),  G(3)  Cost  of  stage  one  and  two,  respectively 

p(2),  P(3)  Aerator  motor  size  of  stage  one  and  two, 

respectively 
§G  Total  cost  of  process 


■ 
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Table  III-l.    (C 

oncinued) 

Symbol 

Explanation 

SP 

Total  horsepower  of  both  aerator  motors 

used  in  the  process 

ST 

Sum  of  holding  times  for  both  stages 

sv 

Sum  of  volumes  for  both  stages 

T(2),  T(3) 

Holding  time  for  stage  one  and  two, 

respectively 

V(2),  V(3) 

Volume  of  stage  one  and  two,  respectively 

X(l),  X(2),  X(3) 

Concentration  of  BOD  going  into  stage  one, 

in  stage  one,  and  in  stage  two,  respectively 

Y 

Y 

• 
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OEL'DMA-B) 
X(2)» A 


Fig.  III-l.   Logic  diagram  for  computer  program  SEARCH. 
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Table   III-2.      Computer  .Program  SEARCH 


■'.CNSS 

•:cnss 

■"CNS  o 
•CNSS 

■<~\SS 
■1CNSS 


.  » ,  ^SEARCH 
I)  ,P(3)  »VC 


COKT   5  MINUTES.   5  PAC 
ASGN  MJB»12 
ASGN  KGC.16 
MODE  GC 
EXEQ  FORTRAN: 
DIMENSIONX(3)»Ct3>»T(3>»G(: 

1  F0RN'AT(7E13«6) 

2  FCPMATl 1HL.13X>3E13.6) 

3  FCRKAT(lHK>26Xt3E13.6) 

4  FORMAT) 1HK»13X»3E13. 6) . 

5  FORMAT! 1H1.7E1 3.6) 

6  R£AD(1»1)AK»DK»CK.Y 
READ(1»11AC1»AC2»AC4»AC5»A2 
READ!! ,1!X( 1) »X(3) .0 
READ(1»1)A»8»0 

WRITEI3.5 )AK,DK.CK»Y 

WRITE<3«1)AC1»AC2»AC4»AC5»A2 

WRITE!3»1)X<  1)  .X(3)  >Q 

WRITE(3»1)A»B»D 

DEL=D*<A-3) 

X (2)=A 
IS  DC2r,\  =  2»3 

C<N)=C(N-l)+( (X(N-l)-X(N) )*(AK*X(N)-DK*(CK+X(N) I )*Y>/< AK*X<N} ) 

T!N>=! <X(N-1 l-X(M) )*Y*(CK+X(N) ) ) / ( AK*X(N>*C <N! ) 

P(N>*Q*(XtN-l!-XtN) ) *2.5958  3E-C5 

V(N)=C*TI N) 

V.'RITE!3»3J  C(N)  »T(N)  »P(N) 

IFfPINi  .LE.20.  1GCTC25 

AC3=?1.57 

Al*i.06 

G0TC27 
2  5  AC3=]0C. 

Al«0.53 
27  GIN)=AC1+AC2*P(N)+AC3*P(N)**A1+AC4*V(N)**A2+AC5*V(N) 

2S  CONTINUE 

SV=V(2!+V(3) 

SP=P(2)+P(3) 

S7=T(2)+T(3) 

SG=G(2)+G(3) 

WRITE(3»2)C(1)»C12> »C(3J 

WRITE (3. 4 1X1 1 ) .XI 2) »X(3) 

WRITE! 3,3 )V( 2) ,V!3) » SV 

WRIT£(3»3)P< 2) »P(3) »SP 

WRITE!  3*3  JT(  2)  »T'(  3)  »ST 

WRIT£(3»3)6< 2) >G(3) »SG 

X!2)=X(2)-DEL 

IF(X(2).LT.3)GCTC6 

S0T0J6 

END 


MONSS 


EXEQ  LINKLCAD 
CALL  SEARCH 
EXEQ  SEARCH. MJ6 


DATA 
.  100OOOE+O0 
■500000E+04 

.  1750CCE+03 
.024000E+03 
.150QCOE+03 
.174990E+03 


.002000E+00 
.149158E+04 

.C01750E+03 
.OlCOOOE+03 
.010000E+03 
.CC1755E+03 


.100000E+03 
.676441E+01 
•451170E+05 
•lOOOOOE+00 
.050000E+00 
.050000E+00 


•50C000E+00 
.009260E+00 


.666667E-I-0C 
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This  work  illustrates  the  application  of  optimization  techniques  to 
water  purification  processes.   Two  different  types  of  water  purification 
processes  are  studied.   A  different  optimization  technique  is  used  to  opti- 
mize each  of  the  processes. 

In  Part  I,  the  Process  Optimization  Program  (POP-II) ,  a  nonlinear  pro- 
gramming technique,  is  used  to  optimize  a  Multi-Effect  Multi-Stage  (MEMS) 
seawater  distillation  process. 

The  results  of  the  optimization  are  in  good  agreement  with  the  results 
obtained  by  applying  a  discrete  version  of  Pontryagin's  Maximum  Principle. 
At  optimum  operating  conditions,  the  one  thousand  gallon  per  hour  MEMS  pro- 
cess optimized,  produces  potable  water  at  a  cost  of  $0.2866  per  one  thousand 
gallons. 

Complex  optimization  problems  with  many  independent  variables  can  be 
rapidly  programmed  for  solution  by  POP-II  because  the  first  partial  deriva- 
tives which  POP-II  requires  are  calculated  by  the  central  differences  numeri- 
cal approximation  technique  rather  than  by  analytic  methods.   However,  some 
manipulation  of  the  input  data  and  several  trials  to  select  the  proper  operat- 
ing mode  may  be  required  before  an  optimum  is  obtained  by  POP. 

In  Part  II,  a  simplified  model  of  a  stage  wise  aerated  lagoon  process 
for  the  partial  purification  of  petroleum  refinery  waste  water  is  derived. 
This  model  was  optimized  using  a  discrete  version  of  Pontryagin's  Maximum 
Principle.   A  method  of  working  with  piecewise  smooth  objective  functions  was 
used  to  obtain  the  optimum  solution. 

The  optimization  study  predicts,  in  agreement  with  theory,  that  tapered 
aeration  is  the  best  mode  of  operation  for  stagewise  operated  aerated  lagoons. 
However,  to  assure  better  quantitative  results,  a  more  comprehensive  aerated 
lagoon  model  is  proposed  for  future  work. 


